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1.  INTRODUCTION:  The  blood-brain  barrier  renders  the  central  nervous  system  a  sanctuary  site  for  disease. 
In  the  era  prior  to  the  development  of  effective  systemic  therapies  this  fact  was  not  of  much  clinical  relevance  as 
most  patients  succumbed  to  visceral  metastatic  disease  before  their  brain  metastases  became  symptomatic. 
However,  as  systemic  breast  cancer  therapies  improve  and  are  able  to  successfully  control  non-Central  Nervous 
System  Disease,  the  brain  is  increasingly  the  first  site  of  relapse.  Both  conventional  chemotherapeutic  agents  and 
targeted  monoclonal  antibodies  do  not  cross  the  blood-brain  barrier  (BBB)  at  concentrations  sufficient  to 
successfully  treat  metastatic  disease.  The  objective  of  the  proposed  work  is  to  use  an  active  form  of  transport  of 
therapeutics  across  the  BBB;  a  transport  mechanism  that  does  not  rely  on  passive  diffusion  or  receptor-mediated 
transcytosis.  We  hypothesize  that  monocytes/macrophages  are  actively  recruited  to  metastases  by  cytokines 
elaborated  by  the  tumors  and  that  these  cells  can  be  utilized  to  transport  therapeutics  directly  to  brain  metastases 
of  breast  cancer.  This  hypothesis  is  supported  by  a  pilot  study  we  published  that  demonstrated  delivery  in  a  mouse 
model  of  fluorescent  microspheres  to  brain  metastases  by  macrophages. (i)  Our  objective  in  this  proposal  is  to 
deliver  a  therapeutic,  doxorubicin,  to  brain  metastases  using  a  nanoshell-double-stranded  DNA  (dsDNA)-drug 
complex  loaded  within  monocytes/macrophages.  Once  present  within  the  metastasis,  we  aim  to  release  the  drug 
by  transcranial  irradiation  at  near-infrared  (NIR)  wavelengths. 

2.  KEYWORDS:  Brain  metastases,  nanoshells,  macrophage,  monocyte,  blood-brain  barrier,  photothermal 
therapy,  light  triggered  release,  Adriamycin,  doxorubicin,  metastatic  breast  cancer,  Trojan  horse 

3.  OVERALL  PROJECT  SUMMARY :  Summarize  the  progress  during  appropriate  reporting  period  (single 
annual  or  comprehensive  final).  This  section  of  the  report  shall  be  in  direct  alignment  with  respect  to  each  task 
outlined  in  the  approved  SOW. 

Task  1 

a.  Current  objectives 

The  original  proposal  suggested  the  immobilization  and  near-IR  triggered  delivery  of  lapatinib  from  gold 
nanoshells.  The  goal  was  to  load  lapatinib  onto  the  nanoshell  surface  through  the  intercalation  of  dsDNA 
sequences  bound  to  the  nanoshell  surface.  The  dsDNA  sequences  would  have  one  strand  that  was  thiolated  and 
covalently  linked  to  the  particle  surface,  while  the  complementary  strand  was  not.  Upon  irradiation,  the 
complementary  sequence  would  be  released  from  the  particle  surface  thereby  releasing  lapatinib  locally.  We 
proposed  to  load  these  nanocomplexes  into  macrophage  /  monocyte  cells  that  could  then  be  used  to  deliver  them 
into  the  brain  metastases. 

In  last  year’s  report,  we  detailed  the  difficulties  we  encountered  with  lapatinib  solubility  and  light- triggered 
release.  As  a  potential  solution,  we  decided  to  focus  on  the  delivery  of  doxorubicin  (DOX).  DOX  is  a 
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chemotherapeutic  agent,  used  in  the  conventional 
treatment  of  cancers  including:  Hodgkin's  lymphoma, 
myeloma,  lung,  ovarian,  gastric,  thyroid,  breast,  sarcoma, 
and  pediatric  cancers.  Its  mechanism  of  action  is  to 
intercalate  between  the  double  strands  of  DNA.  Our 
modified  approach  is  to  intercalate  DOX  into  a  GC  rich 
double-stranded  DNA  (dsDNA)  sequence  and  under  laser 
illumination  at  the  plasmon  resonance  of  the  NS@DNA@DOX  complex,  dehybridize  the  dsDNA  and  release 
DOX  (Scheme  1). 

b.  Results 

In  order  to  deliver  DOX,  we  synthesized  a  nanocarrier  from  which  DOX  could  be  released  by  DNA 
dehybridization.  The  dsDNA  sequences  used  in  this  study  were  engineered  specifically  to  bind  to  the 
chemotherapeutic  DOX:  5’-HS-(CH2)6-CAA  TCA  ATA  GOT  ATC  GTT  CG-3’  and  3’-GTT  AGT  TAT  CGA 
TAG  CAA  GC-5’.  Prior  to  attachment  to  the  nanoshell  (NS)  surface  the  DNA  was  hybridized  by  mixing  equal 
moles  of  the  two  DNA  strands  with  33  mM  NaCl,  heating  the  mixture  to  100  °C  for  4  min  and  slowly  cooling  to 
room  temperature  overnight.  Nanoshells  were  functionalized  with  double-stranded  DNA  (dsDNA)  via  a  Au-S 
bond.  The  dsDNA-NS  mixture  was  incubated  overnight  with  continuous  shaking.  The  solution  was  washed 
several  times  to  remove  the  dsDNA  excess.  The  NS@dsDNA  complex  was  incubated  with  DOX  for  24  h  and 
excess  DOX  was  removed  via  centrifugation. 

Vis-NIR  spectroscopy,  transmission  electron 
microscopy  (TEM),  zeta  potential,  and  dynamic 
light  scattering  (DLS)  measurements  were 
utilized  to  provide  evidence  of  the  formation  of 
the  DOX  nanocarrier.  The  extinction  maximum 
of  bare  NS  slightly  red-shifted  after  dsDNA 
functionalization,  due  to  the  change  in  the 
dielectric  environment  around  the  NS  (Fig.  1A). 

No  changes  in  the  peak  shape  or  width  of  the 
extinction  spectra  was  observed,  indicating  that 
the  NS@dsDNA  nanocomplexes  do  not  aggregate 
during  the  functionalization  process.  TEM  images 
illustrate  the  thin  layer  of  DNA  formed  on  the  nanoshell  surface  (Fig.  IB).  Zeta  potential  measurements  further 
confirm  modification  of  the  NS  surface  (Figure  1C).  The  zeta  potential  of  the  particles  decreased  ~18  mV  after 
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Figure  1.  Formation  and  characterization  of  doxorubicin 
nanocarrier.  (A)  Extinction  spectra  of  NS  (black), 
NS@DNA  (green),  and  NS@DNA@DOX  (red);  (B)  TEM 
images  of  NS  (black)  and  NS@DNA  (green);  (C)  zeta 
potential  of  NS,  NS@DNA,  and  NS@DNA@DOX. 


NS@dsDNA@DOX  NS@ssDNA  ssDNA  +  DOX 

Scheme  1.  Near-IR  light-triggered  release  of 
doxorubicin  (DOX)  from  nanoshells  by 
dehybridization  of  dsDNA. 


DNA  functionalization  as  a  result  of  the  negative  charges  on  the  DNA  phosphate  backbone.  The  zeta  potential 
increased  slightly  after  incubation  with  DOX  due  to  screening  of  the  DNA  charges. 


Cellular  uptake  of  the  NS@DNA  complex  was  investigated  using 
confocal  and  reflectance  (emission  and  excitation  at  640  nm) 
microscopy.  For  the  initial  experiments,  RAW  264.7  a  mouse 
macrophage  cell  line  from  an  Abelson  murine  leukemia  virus-induced 
tumor  is  being.  This  cell  line  has  the  advantage  of  being  immortalized 
and  is  easy  to  propagate  in  cell  culture  in  comparison  to  the  monocytes 
isolated  from  human  blood  and  utilized  previously.  Briefly,  the  RAW 
264.7  macrophage  cells  were  cultured  for  24  h  in  a  p-Slide  8-  well  and 
then  incubated  with  the  NS@dsDNA.  After  20  h,  the  cells  were  rinsed 
three  times  with  PBS  and  fixed  in  a  4%  paraformaldehyde  solution.  The 
cells  were  washed  twice  with  PBS  and  incubated  with  a  Hoechst  33342 
solution  to  stain  the  cell  nucleus.  After  10  min  of  incubation  and 
protected  from  light,  the  cells  were  rinsed  twice  with  PBS.  All  images 
were  collected  on  a  Nikon  Al-Rsi  Confocal  microscope  and  the 
particles  were  imaged  in  the  reflection  mode.  The  3D  maximum 
intensity  projection  image,  corresponding  to  a  slice  from  the  middle  of 
the  cell,  confirms  internalization  of  the  NS@DNA  nanoparticles  (red 
dots,  Fig.  2).  The  details  in  the  orthogonal  YZ  and  XZ  planes  illustrate 
the  location  of  the  nanoparticles  within  the  cell  boundary,  around  but  not  inside  the  nucleus. 


Figure  2.  Confocal  image  of  cellular 
uptake  of  the  nanoshell-DNA 
nanocomplex  in  a  RAW  264.7 
macrophage  cell.  Orthogonal  view 
of  the  3D  maximum  intensity 
projection  image  of  cells  incubated 
overnight  with  nanocarrier.  Nucleus 
is  stained  blue  with  Hoechst  33342. 
Nanoparticles  were  located  by 
reflection  (emission  and  extinction 
at  640  nm).  Scale  bar  =10  pm;  60x 
oil  imersion  objective;  NA  =  1.4. 


In  response  to  the  previous  challenges  we  encountered  with  regard  to  demonstrating  light-triggered  release  of 

lapatinib,  we  decided  to  take  a  step  back  and  investigate  in  greater  detail  the  mechanism  of  light-triggered  release 

in  order  to  better  control  and  tailor  the  release  for  various  molecular  cargos.  In  this  study,  we  used  NS 

functionalized  with  dsDNA;  one  strand  was  thiolated  and  the  complementary  strand  tagged  with  fluorescein 

[5’HS-(CH2)6-TAT  GAT  CTG  TCA  CAG  CTT  GA-3’  and  5’-  ATA  CTA  GAC  AGT  GTC  GAA  CT-Fluor-3’]. 

The  dsDNA  dehybridization  process  was  investigated  by  fluorescence  spectroscopy.  It  had  previously  been 

demonstrated  that  when  a  high  concentration  of  nanoparticles  was  present,  release  of  cargo  using  a  continuous 

wave  (CW)  laser  could  be  achieved  at  temperatures  below  the  DNA  dehybridization  temperature  [1].  Based  on 

this  observation,  a  non-thermal  release  mechanism  by  hot-electrons  was  proposed.  In  this  study,  we  reduced  the 

particle  concentration  to  better  simulate  the  concentration  of  particles  that  is  realistically  achievable  within  a 

tumor.  We  found  by  heating  the  nanoparticle  solution  with  a  temperature-controlled  stage  that  thermally  induced 

dehybridization  of  our  specific  dsDNA  sequence  occurred  at  61.1  ±  2.3  °C  (Fig.  3  A).  Dehybridization  of  the  same 
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sequence  using  CW  light  occurred  at  59.4  ±  5.1  °C  (Figure  3B). 

In  these  experiments,  we  observed  that  release  does  not  occur 
until  the  laser  power  is  high  enough  to  increase  the  bulk 
temperature  above  the  dehybridization  temperature.  This  is 
evidence  that  the  dsDNA  dehybridization  occurs  as  a  result  of 
the  heating  from  the  nanoparticles  as  opposed  to  hot-electrons, 
which  would  be  linearly  dependent  on  the  laser  power.  The  high 
nanoparticle  concentrations  used  in  earlier  experiments  resulted 
in  optical  light  trapping  due  to  multiple  scattering  of  the 
nanoparticles,  allowing  the  local  temperature  to  increase  and 
dehybridize  the  DNA  before  a  significant  rise  in  the  global/bulk 
temperature  was  measured.  Unfortunately,  that  nanoparticle 
concentration  is  not  achievable  in  a  tumor.  Alternatively,  a 
femtosecond  pulsed  wave  laser  enables  significant  temperature 
increases  on  the  nanoparticle  surface  without  heating  the 
solution  (Figure  3C).  By  irradiating  our  nanoparticle  solution 
with  an  ultrafast  laser,  we  were  able  to  demonstrate  release  at 
43. 9±  1.1  °C,~16°C  belowthe  dehybridization  temperature.  By 
further  decreasing  the  spot  size,  we  should  be  able  to  release  the 
DNA  without  any  measureable  temperature  increase. 

Subsequently,  we  examined  the  intracellular  release  of  DOX  by 
both  CW  and  pulsed  lasers.  An  MTT  [3-(4,5-dimethylthiazol-2- 
yl)-2,5-diphenyltetrazolim  bromide]  cell  proliferation  assay  was 
performed  to  investigate  both  the  effects  of  the  nanocarrier  and 
laser  irradiation  on  cell  viability.  RAW  264.7  macrophage  cells 
were  used  for  the  cytotoxicity  studies.  All  experiments  were 
performed  in  triplicate  and  the  percentages  of  viable  cells  were 
calculated  considering  the  absorbance  of  the  control  well  as  a  100%  live  cells.  The  DOX  nanocarrier  was 
incubated  with  cells  overnight,  and  then  irradiated  at  700  mW  for  3  min.  The  cell  proliferation  assay  was 
performed  at  24  h. 

In  all  studies,  an  increase  in  cell  death  was  observed  for  the  DOX  nanocarrier  when  compared  to  cells  with  the 
NS-DNA  nanocomplex  (Fig.  4).  However,  the  same  amount  of  cell  death  was  observed  for  the  cells  with  the 
unilluminated  DOX  nanocarrier  as  for  the  cells  with  the  DOX  nanocarrier  illuminated  with  the  CW  or  pulsed 
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Figure  3.  Comparison  of  the 
dehybridization  temperature  achieved  by 
various  DNA  release  methods;  (A)  thermal 
release,  (B)  continuous  wave  light-triggered 
release,  and  (C)  ultrafast  pulsed  wave  light- 
triggered  release. 


laser.  This  result  suggests  that  the  DOX  may  be  leaking  out  of  the 
complex  before  irradiation.  It  is  theoretically  possible  but 
unlikely  that  the  DOX  nanocarrier  is  toxic  to  the  cells.  We  have 
previously  shown  that  a  similar  nanocarrier  loaded  with  DAPI  is 
not  toxic  to  cells  [2], 

Current  and  future  objectives  include  studies  to  determine  if  DOX 
is  leaking  out  of  the  DNA,  selection  of  an  alternate  drug  that  can 
be  maintained  within  the  dsDNA@NS  complex  until  release  is 
triggered  and  development  of  a  new  host  (i.e.  peptides)  for  the 
DOX  light  triggered  release  process. 


Task  2 

a.  Current  objective 

To  optically  release  the  drug  attached  to  the  surface  of  the  gold  nanoparticle,  it  is  necessary  to  deliver  and  validate 
the  instantaneous  fluence  (the  radiant  energy  per  unit  surface  area)  at  the  location  of  the  tumor.  Given  that,  there 
is  little  data  on  the  amount  and  distribution  of  near  infrared  light  through  the  skull  and  brain  tissue.  To  provide 
an  estimate  of  the  energy  transport  of  near  infrared  light  through  the  skull  and  within  brain  tissue,  a  number  of 
diagnostic  tools  were  developed.  First,  a  nearly  isotropic  optical  dosimetry  probe  was  fabricated  to  obtain  local 
photon  fluence  measurements  with  1 .5  mm  spatial  resolution  and  high  sensitivity.  Calibration  of  the  probe  was 
performed  in  air  and  water-based  phantoms  using  a  calorimeter  (considered  a  Gold  Standard).  Second,  a  voxel- 
based  Monte  Carlo  simulation  code  written  in  CUD  A  for  GTX  Nvidia  (Dr.  Fang,  Harvard)  was  further  developed 
for  different  optical  beam  properties,  such  as  super  Gaussian  beams  with  parallel,  converging,  and  diverging  beam 
profiles,  and  with  the  ability  to  run  on  multiple  Nvidia  cards  (n=3).  Third,  an  empirical  method  to  simulate  photon 
transport  was  developed  to  provide  sufficiently  accurate  solution  at  accelerated  rates.  Monte  Carlo  simulations 
for  photon  transport  are  computational  intensive  and  therefore  are  typically  not  practical  for  clinical  use.  These 
codes  provide  solutions  in  real-time  during  the  developmental  stage  of  the  treatment  plan.  These  tools  were  used 
to: 
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Figure  4.  Intracellular  light-triggered 
release  of  doxorubicin.  MTT  plate-based 
cell  viability  assay  of  RAW  264.7  mouse 
macrophage-like  cells  alone  (blue)  and  cells 
incubated  with  NS@DNA  (green)  and 
NS@DNA@DOX  (red)  for  24  h.  Cells  were 
irradiated  with  a  CW  or  ultrafast  pulsed 
laser  (700  mW,  3  min).  Cell  viability  was 
evaluated  24  h  later. 


2a.  Estimate  of  the  attenuation  of  near  infrared  irradiation  by  the  skull  using  the  fluence  probe;  and 
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2b.  Estimate  of  the  attenuation  of  near  infrared  irradiation  by  the  skull  and  brain  tissue  phantom  using  the 
fluence  probe. 


b.  Results 

i.  Calibration  of  the  Optical  Dosimetry  Probe  (ODP)  [1,3]:  The  objective  was  to  fabricate  and  calibrate  the 
optical  dosimetry  probe  (ODP)  to  a  calorimeter  (Nova  II;  Ophir).  The  ODP  consists  of  a  1.5  mm  diameter  sphere 
(epoxy  and  titanium  dioxide)  molded  onto  the  end  of  a  200pm  optical  fiber.  The  output  of  the  fiber  was  coupled 
to  a  photodiode  (ET-2030,  Electro  Optics  Technology  Inc),  the  output  of  which  was  recorded  by  a  digital 
oscilloscope  (Tektronix  TDS3052B).  The  digitized  signal  was  uploaded  to  a  PC,  where  the  signal  was  passed 
through  a  high-pass  butterworth  filter  to  reduce  background  noise  and  the  area  underneath  the  pulse  calculated  to 
determine  the  photon  fluence.  To  calibrate  the  ODP,  the  photon  fluence  in  air  was  measured  by  both  the  ODP  and 
calorimeter  using  a  uniform  laser  beam  formed  by  passing  the  output  of  the  laser/OPO  system  (Nd:YAG  class  4 
laser/  MagicPRISM  OPO,  Opotek  Inc.)  (tuned  to  800nm)  through  an  integrating  sphere  (Melles  Griot).  The  results 
are  displayed  in  Figure  1.  To  determine  the  fluence  in  water  (or  tissue),  the  change  in  the  index  of  refraction 
causes  the  reflective  properties  between  the  probe  media  and  water  to  change,  thus  requiring  an  additional 
correction  factor.  The  probe  was  placed  within  a  macrocuvette,  a  plexiglass  container  which  was  approximately 
5.0  x  5.0  x  5.0  cm3  (5mm  thick  walls).  The  ODP  signal  and  calorimeter  measurements  were  acquired  in  air  as 
before.  The  macrocuvette  was  then  filled  with  water,  and  the  ODP  was  positioned  along  the  central  axis  of  the 
uniform  beam.  From  this  data,  a  conversion  factor  was  applied  to  the  ODP  signal  to  convert  the  calibration  curves 
in  air  (Figure  1)  to  tissue.  This  conversion  factor  also  took  into  account  the  variation  in  the  beam  width  when 
traversing  from  plexiglass-to-air  versus  plexiglass-to-water  and  photon  attenuation  for  optical  path  lengths  within 
water,  both  of  which  were  relative  small  (~3%,  the  average  error  in  the  laser  output  power). 
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Figure  1 :  (Left)  The  output  of  the  ODP  in  mVolts  plotted  as  a  function  of  laser  fluence  measured  by  a 
calorimeter.  (Right)  The  angular  dependency  (azimuthal  and  rotation)  of  the  ODP  to  the  direction  of  laser 
light  was  calculated  as  a  function  of  the  absolute  value  of  the  residuals  and  how  it  compares  to  nylon 
probes. 


ii.  Phantoms  Reproducing  the  Optical  Properties  of  the  Head:  The  objective  was  to  develop  a  formula  that 
could  be  used  to  create  water-based  (and  agar-based)  phantoms  with  the  desired  optical  properties  (pa  and  ps’) 
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for  the  skull,  white  matter  and  gray  matter  (at  800nm,  the  wavelength  of  stimulated  drug  release).  India  ink  was 
used  as  the  absorber  (with  some  but  little  scattering),  while  intralipid  was  used  as  a  scatterer  (with  some  but  little 
absorption).  To  determine  the  absorption  coefficient,  various  concentrations  of  india  ink  in  water  were  placed 
within  the  photospectrometer,  from  which  a  calibration  curve  plotted:  15.07  x  %  (india  ink)  +  0.0049  (^=0.99) 
(data  not  shown).  The  reduced  scattering  coefficient  was  determined  by  placing  different  concentrations  (from  1 
to  7%)  of  intralipid  in  water,  exposing  the  solution  to  a  uniform  (1cm  diameter)  laser  beam  using  an  integrating 
sphere  as  described  in  the  previous  section.  The  photon  fluence  was  measured  using  the  ODP  at  various  path 
lengths  (or  depths)  along  the  central  axis  of  the  beam.  To  determine  the  reduced  scattering  coefficient  in  each 
sample,  Monte  Carlo  simulations  were  performed  for  various  reduced  scattering  coefficients  until  a  minimum  in 
the  chisquare  was  reached.  The  data 
was  plotted  and  fit  to  a  line:  11.30  x 
(%  intralipid)  -  2.17  (r2=0.99)  (data 
not  shown).  From  these 
measurements,  the  white  and  gray 
matter  phantoms  were  fabricated 
based  the  optical  properties 
measured  by  Van  der  Zee  et  al  [6],  and  the  skull  phantom  optical  properties  were  based  on  measurements  by 
Firbank  et  al  [7],  see  Table  I. 

iii.  Simulating  Photon  Transport  in  Phantoms/Tissue  [  1  -3] :  The  objective  was  to  develop  and  validate  a  voxel- 
based  Monte  Carlo  code  of  photon  transport  in  the  head  and  brain.  A  voxel-based  3-D  Monte  Carlo  (MC)  code 
(Fang  et  al  [8])  developed  in  CUDA  for  NVidia  cards  was  modified  to  allow  the  user  to  create  laser  sources  with 
various  aperture  profiles,  such  as  uniform  and  super-Gaussian,  and  directionality,  such  as  diffusive,  converging, 
and  diverging  beams.  Supporting  code  (Matlab)  was  developed  to  allow  the  user  to  design  their  own  objects, 
which  in  this  study  was  based  on  phantom  designs  and  medical  images  (CT). 

A  critical  advantage  of  a  voxel-based  code  is  not  only  to  include  tissue  heterogeneities  but  to  integrate  medically 
acquired  images  into  the  simulations.  To  date,  the  MC  code  has  not  been  validated  with  measured  data;  therefore, 
a  number  of  experiments  were  performed.  First,  the  photon  fluence  distribution  was  measured  using  the  ODP  in 
phantoms  representing  white  matter,  gray  matter,  and  skull  bone.  The  ODP  probe  was  secured  to  a  3-D  stage 
micrometer  to  obtain  measurements  along  the  central  axis  of  the  beam  at  various  depths  in  the  phantom  (axial 
profile)  and  perpendicular  to  the  central  axis  (lateral  profile).  These  measurements  were  taken  using  a  uniform 


Table  I:  Optical  properties  of  the  skull  bone,  white  matter,  and  gray 
matter. 
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550 
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0.24 

184 

0.9 

18.4 

1.82 

0.0156 
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and  a  super  Gaussian  laser  beam.  Next,  Monte  Carlo  simulations  were  performed  and  compared  to  measurement. 

In  Figure  2,  the  normalized  fluence  along  the  central  axis  of  the  flat  beam  was  plotted  for 

both  white  and  gray  and  for  ODP  and  MC  data.  For  depths  up  to  2  cm,  the  agreement  is  within  3%.  Similar 


experiments  were  performed  for  the  super  Gaussian  beam,  where  for  depths  up  to  3  cm,  the  MC  simulated  results 
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Figure  2.  The  fraction  of  the  input  fluence  as  measured  by  the  ODP  at  different  depths  in  white  and  gray  matter 
phantoms  for  a  uniform  laser  beam  (1.0cm  in  diameter)  (left)  and  a  super  Gaussian  beam  profile  (right).  Included 
in  these  plots  are  the  Monte  Carlo  simulated  results,  which  agree  to  within  3.0%  of  the  measured  values. 


matched  the  ODP  measurements.  The  comparison  of  the  lateral  profiles  for  the  super  Gaussian  beam  is  shown  in 
Figure  3,  and  demonstrates  the  excellent  agreement  between  ODP  and  MC  results. 


Figure  3.  Plotted  are  the  relative  photon  fluence  as  measured  by  the  ODP  and  as  simulated  by  the  MC  code  for  a  super  Gaussian 
laser  beam  in  air  (top  left),  white  matter  (top-right),  gray  matter  (bottom -I eft),  and  skull  (bottom-right). 
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iv.  Sheep  Skull  and  Brain  Phantom:  In  combination  with  the  phantom  study  from  the  previous  section,  an 
estimate  of  the  attenuation  of  near  infrared  irradiation  by  the  skull  and  brain  tissue  phantom  was  made.  The  laser 
beam  emanating  from  the  optical  fiber  bundle  was  measured  in  air  using  the  ODP  and  used  to  determine  the  beam 
profile,  e.g.,  super  Gaussian  function,  and  laser  fluence.  This  laser  beam  was  positioned  as  close  to  normal  to  the 
goat  skull  as  possible.  The  thickness  of  the  bone  where  the  laser  beam  was  exposed  was  measured  to  be  1.1cm 
using  calipers  and  validated  by  CT  scans.  To  measure  the  laser  beam  profile  at  the  inner  surface  of  the  skull,  the 
dosimetry  probe  was  axially  positioned  within  4  mm  of  this  inner  surface,  and  using  a  3-D  micrometer,  the  ODP 
was  translated  perpendicular  to  the  laser  beam  to  measure  the  fluence  profile  and  along  the  axis  of  the  beam 
profile  to  determine  the  peak  beam  intensity  as  a  function  of  the  distance  from  the  inner  surface  (see  Figure  4). 
After  applying  the  calibration  factors  (Section:  Calibration  of  ODP),  the  percent  of  the  laser  fluence  transmitted 
through  the  skull  was  measured  to  4.41%  (SA2.A).  MC  simulated  data  were  consistent  with  these  ODP 
measurements  (see  Figure  4). 
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Figure  4.  (Top  row)  Displayed  are  the  photos  of  the  goat  skull  and  the  ODP  used  in  the  phantom  experiments.  The  picture  on 
the  right  has  the  laser  beam  exposing  the  frontal  bone  (left  to  right)  and  the  ODP  embedded  in  the  brain  tissue  phantom.  (Bottom 
row)  Plotted  are  the  lateral  profile  near  the  inner  surface  of  the  skull  and  the  axial  profile  along  the  central  axis  of  the  laser  beam 
as  measured  by  the  ODP  probe  and  simulated  using  the  MC  code. 


With  this  same  setup,  the  skull  cavity  was  filled  with  a  4%  intralipid  solution.  The  ODP  was  used  to 
measure  the  attenuation  of  the  beam  (Figure  4.c):  the  transmitted  fluence  was  measured  to  be  7.25%  per  cm, 
which  is  consistent  with  Figure  2.  If  a  7.5%  of  intralipid  (white  matter)  was  used,  the  measured  transmitted  fluence 
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was  3.1%  per  cm,  and  gray  matter  was  measured  to  be  0.3 1%  per  cm  (see  Figure  2).  To  relate  these  measurements 
to  a  human  patient,  average  values  for  the  thickness  of  the  cortical  gray  matter  and  skull  bone  is  determined. 
Cortical  gray  matter  thicknesses  can  range  from  2.0  to  5.0  mm,  averaging  around  4.7-4.85  mm  [9,  10].  Skull  bone 
thickness  can  range  from  5.0  to  1 1.0  mm  in  an  adult  human  head,  averaging  around  7.8  mm  [11],  Based  on  the 
above  measurements,  the  average  fluence  transmitted  through  8.0  mm  of  the  skull  bone  would  be  10.3% 
and  through  4.8mm  of  gray  matter  would  be  9.1%,  or  0.93%  on  average  through  both  tissues.  Therefore, 
a  tumor  4.5cm  deep  within  the  brain  would  receive  approximately  8.5xl0~9  of  the  laser  fluence  (SA2.B1 
Depending  on  the  location  of  the  tumor,  the  positioning  of  the  optical  sources  could  take  into  account  the 
heterogeneity  of  the  skull  and  gray  matter.  Optimally,  a  region  of  the  skull  and  gray  matter  that  that  is  5mm 
(24.3%  transmitted)  and  2.5mm  (20.0%  transmitted)  thick  would  increase  the  fluence  to  1.9x1  O'8,  an  improvement 
of  a  factor  of  4.4.  This  suggests  that  a  multiple  number  of  sources  optimally  positioned  for  treatment  based  on 
the  patient’s  anatomy  may  be  beneficial. 

v.  Empirical  Model  of  Photon  Transport  [2,3]:  An  empirical  method  was  developed  to  simulate  photon 
transport  in  tissue  based  on  the  acquired  medical  images  of  an  individual  patient  (CT  and  MRI)  with  sufficient 
accuracy  and  at  accelerated  rates  for  clinical  use.  This  approach  is  based  on  the  assumption  that  the  photon  fluence 
in  a  voxel  within  a  particular  layer  (or  slice  of 
an  image)  is  the  weighted  sum  of  the  fluence 
from  its  neighboring  voxels  in  the  previous 
layer.  These  weights  determine  the  fraction  of 
photons  within  a  photon  packet  that  is 
transported  from  one  voxel  to  its  neighbors 
relative  to  its  average  (or  weighted)  direction 
of  propagation.  It  depends  not  only  on  a 
voxel-based  geometry  of  the  image  but  the 
voxel  size  and  optical  properties  of  the 
medium.  It  has  been  shown  to  faithfully 

reproduce  the  local  shadowing  effects  from 

,  .  ,  ,  ,  .  .  ,  ,  ,  ,  ,  ,  •  ,  Figure  5.  Displayed  is  a  diagram  depicting  the  local  photon  transport 

high  absorbers  (e.g.,  blood  vessels)  and  high  of  the  empirical  model  for  a  voxel  within  tissue.  Note,  this  takes  into 

^  ^  \  i  account  the  general  direction  of  the  photon  packet  and  the  scattering 

contrast  media  (Gold  nanoparticles),  and  properties  of  the  tissue. 

general  attenuation  in  homogeneous  tissue 

types.  To  determine  the  weights,  MC  simulations  were  performed  in  known  media,  such  as  white  matter,  gray 

matter,  skull  bone,  and  an  Astrocytoma  (tumor).  Using  these  fluence  profiles,  the  weights  were  modified  until  the 

photon  fluence  distribution  from  the  empirical  model  matched  (minimum  in  Chisquare)  MC  simulated  results  for 

different  voxel  sizes.  A  comparison  between  the  photon  fluence  within  each  individual  tissue  type  and  for  a 
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software  phantom  of  the  skull  and  brain  are  shown  in  Figure  6.  There  remains 
(over  estimate)  of  the  photon  fluence  which 
requires  further  investigation,  but  the 
agreement  is  very  good  to  develop  and 
optimize  delivery  systems  based  on  an 
individual  patient’s  anatomy.  When 
implemented  on  CUDA,  the  time  it  takes  to 
simulate  a  human  head  with  a  broad  beam 
source  (e.g.,  super  Gaussian  beam)  was  on  the 
order  of  seconds  to  low  minutes  where  the  MC 
would  run  in  excess  of  24  hours. 


a  relative  small  systematic  effect 
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Future.  We  are  currently  applying  our 
methods  to  mouse  and  dog  models,  and  human 
scans.  Multiple  sources  (lasers,  semiconductor 
arrays)  coupled  with  interventional  methods 
are  being  investigated  to  improve  NIR  delivery 
to  the  brain.  However,  other  methods  will  also 
be  explored,  such  as  radiation  sensitive  films  to  release  these  drugs. 
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Figure  6.  Comparison  of  the  Empirical  algorithm  with  MC 
simulations  in  homogenous  (a, b)  and  heterogeneous  (c,d) 
phantoms  that  resemble  white  matter,  gray  matter  and 
astrocytoma.  The  plots  in  panels  (a),  (b),  and(d)  display  the 
exponential  decrease  in  fluence  along  the  central  axis  of  a 
Gaussian  beam. 


Tasks  1  d  &e,  2c,  and  Task  3:  These  goals  were  unable  to  be  met  because  of  the  issues  with  drug  solubility 
(laptinib),  light-triggered  release,  and  probable  leakage  of  doxorubicin  from  the  NS@DNA@DOX  nanocarriers. 
They  all  required  the  delivery  of  NS@DNA  with  stably  loaded  drug  (DOX),  which  could  be  released  by  laser 
irradiation. 


4.  KEY  RESEARCH  ACCOMPLISHMENTS:  Bulleted  list  of  key  research  accomplishments  emanating  from 
this  research.  Project  milestones,  such  as  simply  completing  proposed  experiments,  are  not  acceptable  as  key 
research  accomplishments.  Key  research  accomplishments  are  those  that  have  contributed  to  the  major  goals 
and  objectives  and  that  have  potential  impact  on  the  research  field.  If  there  is  nothing  to  report,  simply  state 
"Nothing  to  report. " 

•  Nanocarriers 

Earlier  experiments  by  the  Rice  Nanophotonics  Group  suggest  that  DNA  dehybidization  and  cargo  release 

occurred  by  a  non-thermal  release  mechanism  involving  hot-electrons  [1],  The  collective  oscillation  of 

surface  conduction  electrons  in  metallic  nanoparticles  induced  by  light  at  a  resonant  frequency  is  known  as  a 

localized  surface  plasmon.  Upon  excitation,  surface  plasmons  can  decay  radiatively  by  re-emitting  photons 
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of  light  or  non-radiatively  by  the  generation  of  local  heat  and/or  energetic  electrons.  These  energetic  electrons 
are  known  as  “hot-electrons”  and  the  Rice  group  had  hypothesized  that  the  two  DNA  strands  on  the  NS 
surface,  which  are  negatively  charged,  were  dehybridized  by  electrostatic  interaction  with  the  hot-electrons. 
It  is  now  clear  from  the  experiments  conducted  as  a  part  of  this  research  that  at  particle  concentrations 
achievable  within  a  tumor,  dehybridization  of  the  DNA  double  helix  is  a  thermal  phenomenon.  This  resolves 
a  fundamental  question  as  to  how  the  NS@DNA  nanocarriers  function  to  release  cargo.  It  has  also  led  to  the 
testing  of  the  femtosecond  pulsed  laser  as  a  work-around,  which  will  enable  sufficient  heating  at  the 
nanoparticle  surface  to  release  cargo  but  insufficient  bulk  heating  to  cause  protein  denaturation  and  subsequent 
tissue  damage. 

•  Delivery  of  infrared  irradiation  to  brain 

As  pointed  out  earlier,  prior  to  our  work  there  was  little  data  on  the  amount  and  distribution  of  near  infrared 
light  able  to  be  delivered  through  the  skull  and  brain  tissue.  We  have  been  able  to  model  this  using  both  the 
sheep  skull  and  phantom  as  well  as  the  Monte  Carlo  simulations.  This  has  led  us  to  proposing  that  multiple 
simultaneous  sources  of  IR  irradiation  with  or  without  intracranial  intervention  methods  is  likely  to  be 
required  to  deliver  a  sufficient  fluence  of  irradiation  to  effect  drug  release. 

5.  CONCLUSION:  Summarize  the  importance  and/or  implications  with  respect  to  medical  and  /or  military 
significance  of  the  completed  research  including  distinctive  contributions,  innovations,  or  changes  in  practice  or 
behavior  that  has  come  about  as  a  result  of  the  project.  A  brief  description  of future  plans  to  accomplish  the  goals 
and  objectives  shall  also  be  included. 

The  delivery  of  effective  therapeutics  across  the  blood-brain  barrier  remains  as  critical  a  clinical  challenge  as  it 
was  when  we  first  proposed  our  study.  We  remain  convinced  that  the  active  delivery  by  macrophages  of 
nanobound  therapeutics  to  intracerebral  malignancies  is  a  promising  strategy.  Our  studies  have  provided 
fundamental  information  on  a  number  of  aspects  of  the  delivery  of  therapeutics  by  NS@DNA.  1.  It  is  of  key 
importance  that  the  therapeutic  is  soluble  and  stable  in  aqueous  solution;  2.  Dehybridization  of  DNA  duplexes 
on  the  nanoshell  surface  is  a  thermal  phenomenon.  3.  Although  doxorubicin  is  a  DNA  intercalating  agent,  it  may 
distort  the  double  helix  such  that  it  is  not  a  secure  cargo.  4.  We  have  estimated,  to  the  best  of  our  knowledge  for 
the  first  time,  the  amount  and  distribution  of  near  infrared  light  able  to  be  delivered  through  the  skull  and  brain 
tissue.  While  we  have  demonstrated  substantial  attenuation  of  the  IR  irradiation  on  its  pathway  through  the  skull 
and  brain  phantom,  our  results  have  prompted  us  to  consider  utilizing  multiple  NIR  sources.  Although  this  marks 
the  end  of  our  DoD  funding,  we  hope  to  continue  to  pursue  this  novel  delivery  method.  Should  we  secure  future 
funding,  we  will  determine  if  and  to  what  extent  DOX  is  leaking  from  the  DNA  duplex;  presuming  that  DNA 
duplex  distortion  is  responsible  for  the  leakage,  we  will  search  for  an  alternate  drug  that  can  be  maintained  within 
the  dsDNA@NS  complex  until  release  is  triggered;  and  we  will  begin  development  of  a  new  host  (i.e.  peptides) 
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for  the  DOX  light  triggered  release  process.  With  regard  to  delivering  sufficient  IR  irradiation  to  the  brain,  the 
following  are  strategies  that  we  hope  to  pursue: 

(1)  Place  multiple  sources  to  avoid  thicker  regions  of  gray  matter  and  overall  skull  density. 

(2)  Increase  the  number  of  sources  around  the  skull  in  order  to  get  the  fluence  to  plateau  (instead  of 

exponentially  decay)  at  larger  depths  in  the  brain. 

(3)  Examine  the  possibility  of  directional  scattering  in  neuronal  fibers  as  a  means  to  increase  NIR  depth 
penetration. 

(4)  Given  the  diffusive  nature  of  light,  it  turns  out  that  the  beam  surface  area  can  be  increased  to  reduce 

surface  heating  near  the  surface  of  the  brain  without  too  much  loss  in  the  fluence  at  depth. 

We  will  also  explore  other  methods  such  as  radiation  sensitive  films  from  which  to  release  drugs. 
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a.  List  all  manuscripts  submitted  for  publication  during  the  period  covered  by  this  report  resulting  from  this 
project.  Include  those  in  the  categories  of  lay  press,  peer-reviewed  scientific  journals,  invited  articles,  and 
abstracts.  Each  entry  shall  include  the  author(s),  article  title,  journal  name,  book  title,  editors(s),  publisher, 
volume  number,  page  number (s),  date,  DOI,  PMID,  and/or  ISBN.  Journal  publications 
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2.  Peer-Reviewed  Scientific  Journals: 

[1]  Verleker  AP,  Fang  Q,  Choi  M-R,  Clare  S,  Stantz  KM.  An  empirical  approach  to  estimate  near-infra-red  photon 
propagation  and  optically  induced  drug  release  in  brain  tissues.  Proc  SPIE  BiOS  9308-29:  1-8,  2015  (doi: 
10.1117/12.2079991) 

[2]  Verleker  AP,  Fang  Q,  Choi  M-R,  Clare  S,  Stantz  KM.  An  Optical  Therapeutic  Protocol  to  treat  brain 
metastasis  by  mapping  NIR  activated  drug  release:  a  Pilot  Study.  Proc  IEEE  MIC  (M19-96):  107,  2014. 

[3]  Verleker  AP,  Fang  Q,  Choi  M-R,  Clare  S,  Stantz  KM.  Dosimetric  Validation  of  a  Monte  Carlo  based  Optical 
Planning  Tool,  soon  to  be  submitted  to  Applied  Optics 
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8.  REPORTABLE  OUTCOMES: 
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development  of  cell  lines,  tissue  or  serum  repositories,  funding  applied for  based  on  work  supported  by  this  award, 
and  employment  or  research  opportunities  applied  for  and/or  received  based  on  experience/training  supported 
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This  project  has  directly  contributed  to  the  research  thesis  of  Akshay  Verleker  at  Purdue  University,  who  will  be 
graduating  with  a  doctorate  in  2016.  It  provided  him  the  opportunity  to  submit  and  obtain  a  Purdue  Research 
Foundation  Grant  (2014-2015),  a  Center  for  Cancer  Research  Travel  Award  (2014),  and  the  Graduate  Student 
Travel  (Compton)  Award  from  the  College  of  Health  and  Human  Sciences  (HHS)  (2015).  These  latter  awards 
provided  Akshay  the  opportunity  to  present  this  research  at  the  IEEE  Medical  Imaging  (NSS/MIC  2014)  and 
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photon  transport  and  imaging  sciences. 
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An  Optical  Therapeutic  Protocol  to  treat  brain  metastasis 
by  mapping  NIR  activated  drug  release:  A  Pilot  Study 

Akshay  Prabhu  Verleker,  Qianqian  Fang,  Mi-Ran  Choi,  Susan  Clare  and  Keith  M.  Stantz 


Abstract-Treatment  of  Central  Nervous  System  (CNS) 
metastasis  poses  a  critical  clinical  challenge  due  to  limitations  in 
drug  uptake  across  the  blood  brain  barrier  and  blood- 
cerebrospinal  fluid  barrier.  Recent  research  has  shown  the 
efficacy  of  using  macrophages  as  drug  carriers  to  target 
metastatic  sites  in  the  brain,  which  can  then  be  activated  by 
illuminating  with  near  infrared  radiation.  The  goal  of  this 
research  is  to  develop  an  optically  targeted  therapeutic  treatment 
of  metastasis,  specifically  for  metastatic  breast  cancer  of  the 
brain.  As  a  first  step  towards  accomplishing  this,  we  developed  a 
3D  Monte  Carlo  photon  transport  code  capable  of  simulating 
phantoms  with  optical  properties  of  brain  and  tumor;  and 
validated  the  Monte  Carlo  simulated  photon  fluence  within  brain 
phantoms  using  an  optical  dosimetry  probe.  The  phantom  studies 
showed  good  correlation  (correlation  coefficient  R=0.977) 
between  the  probe  measurements  and  the  Monte  Carlo 
simulation  in  a  white  matter  phantom  (reduced  scattering 
coefficient  ps^S.lSmm"1,  absorption  coefficient  pa=0.005  mm'1). 
Our  future  steps  will  be  to  implement  the  Monte  Carlo  to  map 
out  photon  energy  distribution  in  the  brain,  and  subsequent  drug 
release,  by  segmenting  &  translating  head  CT  image  volumes  to 
corresponding  optical  properties  of  brain  tissues.  To  access  the 
therapeutic  response,  changes  in  the  vascular  physiology  of  the 
brain  due  to  Her2  inhibition  will  be  measured  using  dynamic 
contrast-enhanced  imaging  (e.g.  DCE-CT),  and  with  Monte  Carlo 
based  optical  fluence  maps.  An  optical  treatment  plan,  using  fast 
Monte  Carlo  software,  optimized  with  CT  segmented  image 
volumes,  would  significantly  reduce  the  treatment  time  and  allow 
targeted  drug  activation  while  sparing  healthy  tissues. 

I.  Introduction 

reatment  of  Central  Nervous  System  (CNS)  metastasis 
poses  a  critical  clinical  challenge  due  to  limitations  in  drug 
uptake  in  the  brain  (across  the  blood  brain  barrier  and  blood- 
cerebrospinal  fluid  barrier)  and  adverse  neurotoxic  effects  of 
mainstay  therapies  such  as  whole  brain  radiation  therapy 
(WBRT)  and  stereotactic  radiosurgery.1  With  non-CNS 
metastasis  being  treated  successfully  with  receptor  targeted 
drug  delivery,  aided  by  improved  imaging  and  localized 
radiotherapy  techniques,  the  brain  remains  a  sanctuary  for 
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metastasis  (e.g.  triple  negative  and  HER2+  breast  cancer 
metastasis  with  median  survival  rates  of  2  to  16  months).1  This 
study  is  part  of  a  project  which  demonstrated  a  macrophage 
based  “Trojan  Horse”  delivery  of  fluorescent  molecules  into 
brain  metastasis2  and  subsequent  release  of  the  fluorescent 
molecules  on  activation  in-vivo  by  NIR  radiation.3  The  larger 
goal  of  this  collaborative  project  is  to  deliver  lapatinib-gold 
nanocomplex  laden  macrophages  to  brain  metastasis  and 
design  a  therapeutic  protocol  to  optically  stimulate  drug 
release  in  target  tissues  in  the  brain.  In  the  current  project,  the 
purpose  is  to  design  an  optical  simulation  protocol  (using 
Monte  Carlo  simulations)  to  effectively  predict  the  photon 
distribution  and  subsequent  drug  activation  in  the  brain  tissue, 
and  validate  it  by  using  an  optical  dosimetry  probe.  With  the 
rate/quantity  of  drug  release  being  directly  dependent  on  the 
optical  power  delivered,  the  Monte  Carlo  simulation  protocol, 
with  the  ability  to  predict  photon  distribution  in  heterogeneous 
brain  tissue,  can  be  used  to  generate  an  equivalent  drug  release 
map.  The  optical  properties  of  the  brain  can  be  mapped  by 
converting  the  CT  image  densities  of  the  brain  into  the 
corresponding  optical  properties  of  the  tissue  (e.g.  white/gray 
matter,  CSF,  etc.),  which  would  serve  as  an  input  to  the  Monte 
Carlo  routine.  The  optical  drug  release  map  generated  by 
Monte  Carlo  can  be  co-registered  with  CT  images  to  improve 
the  therapeutic  efficacy  of  the  protocol. 

In  order  to  predict  the  photon  energy  distribution  in 
heterogeneous  brain  tissue,  we  used  a  GPU  based  3D  Monte 
Carlo  simulation  code4,  which  was  modified  to  generate  user 
defined  beam  profiles  to  emulate  the  NIR  (800  nm)  laser 
source  used  in  our  lab.  Monte  Carlo  simulation  has  been  used 
as  the  gold  standard  for  photon  propagation  studies  (in  3D 
heterogeneous  media),  with  the  GPU  version  having 
acceleration  speeds  100  to  300  times  faster  than  corresponding 
CPU  versions.4  In  order  to  validate  and  quantify  the  Monte 
Carlo  generated  fluence  in  brain  tissue  phantoms,  we  used  an 
optical  dosimetry  probe  to  measure  the  photon  fluence  in  low 
optical  fluence  regions.  An  optimized  treatment  plan  using  a 
fast  Monte  Carlo  software  to  predict  in  vivo  drug  release, 
would  significantly  reduce  the  treatment  time  and  allow  the 
targeted  drug  activation  while  sparing  healthy  tissues. 

In  addition  to  designing  the  Optical  therapeutic  protocol,  we 
also  looked  into  the  feasibility  of  using  Dynamic  Contrast 
Enhanced  CT  (DCE-CT)  imaging  to  monitor  and  quantify  the 
efficacy  of  treatment.  The  vascular  structure  of  a  tumor  is 
complex  and  heterogeneous,  and  is  seen  to  vary  between 
individual  tumors  and  during  different  phases  of  tumor 
growth.5  This  heterogeneity  in  tumor  vasculature, 
characterized  by  hypoxia,  anoxia  and  hypoglycemia,  leads  to 
differing  response  to  treatment  due  to  formation  of  barriers  to 
transport  and  distribution  of  drug  molecules.5’9  Following 
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treatment  with  anti-angiogenic  or  cytotoxic  drugs,  the  tumor 
vasculature  undergoes  complex  physiological  and  biological 
changes  some  of  which  are  adaptive  in  nature.5’9  DCE-CT  has 
been  successfully  used  to  simultaneously  quantify 
physiological  parameters  such  as  perfusion  and  vascular 
plasma  volume  in  tumors  after  treatment  with  therapeutic 
drugs.10'12  Past  research  has  proved  the  feasibility  of  using 
DCE-CT  in  measuring  the  physiological  response  of  breast 
tumors  to  anti-angiogenic  therapies  using  the  two- 
compartmental  model.11’12  In  the  present  study  we  used  the 
above  model  to  obtain  parametric  maps  of  perfusion, 
permeability,  fractional  plasma  volume  and  fractional 
interstitial  volume  in  glioma  tumors  grown  in  mice;  thus 
demonstrating  the  feasibility  of  DCE-CT  imaging  to  monitor 
brain  tumor  response  to  treatment. 

II.  Methods 

Step  1:  The  optical  dosimetry  probe  (1.5mm  diameter)  was 
calibrated  using  a  calorimeter  (Gold  standard)  in  low  fluence 
regions.  The  probe  measurements  were  used  to  validate  Monte 
Carlo  generated  photon  fluence  in  brain  phantoms  resembling 
white  and  gray  matter,  made  with  specific  concentrations  of 
India  ink  (absorber)  and  Intralipid  (scatter),13'14  as  measured 
by  photo-spectrometer  (Gold  standard),  consistent  with  brain 
tissue  and  tumors  (values)  as  shown  in  Fig  1.  The  source  of 
light  used  was  an  integrating  sphere  connected  to  a  pulsed 
NIR  laser  with  wavelength  800  nm,  frequency  200MHz,  and 
energy  lOmJ/pulse.  The  broad  beam  NIR  output  power 
emanating  from  the  integrating  sphere  was  0.2083mW/mm2. 

Step  2:  In  order  to  apply  the  Monte  Carlo  (MC)  simulation 
protocol  to  the  brain,  we  designed  a  segmentation  software  to 
determine  tissue  boundaries  in  the  brain  using  Computed 
tomography  images.  The  CT  image  densities  in  Hounsfield 
units  (HU)  of  the  brain  were  translated  into  corresponding 
optical  properties  (absorption  coefficient,  scattering 
coefficient,  and  anisotropy  factor)  using  lookup  tables  (Fig  2). 
The  lookup  tables  were  prepared  from  earlier  studies  which 
derived  the  Hounsfield  units  and  optical  properties  of  brain 
tissues  such  as  white/gray  matter.13'15  The  segmentation 
software  generated  an  optical  map  of  the  brain,  which  was 
given  as  input  to  the  Monte  Carlo.  The  photon  energy  map 
generated  by  the  Monte  Carlo  can  be  co-registered  with  the 
CT  image  to  visualize  the  efficacy  of  drug  release. 

Step  3:  “Optical  treatment  protocol”,  voxel-based  optical 
properties  of  a  human  head  (obtained  from  CT  scans)  were 
input  into  MC  code  to  simulate  3D  fluence  maps  from 
multiple  (optically  coupled)  light  sources.  This  study  also 
looked  into  the  feasibility  of  monitoring  effectiveness  of 
targeted  drug  release  by  monitoring  the  physiology  at  target 
sites  and  surrounding  tissues.  An  Orthotopic  model  of 
glioblastoma  tumors  were  grown  in  mice  and  parametric 
images  of  perfusion,  permeability,  fractional  plasma  volume 
and  fractional  interstitial  volume  were  derived  using  the  two 
compartmental  model.12 
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Fig.  1.  Validation  of  Monte  Carlo  simulated  photon  fluence  (l/mm2)  with 
probe  measurements  (peak  area  in  mV)  in  a  white  matter  phantom.  Log 
normalized  fluence  plotted  v/s  depth  in  a  white  matter  phantom.  Correlation 
R=l. 


III.  Observations  &  Results 


The  Monte  Carlo  generated  photon  fluence  was  validated 
by  probe  measurements  in  optical  brain  phantoms  (e.g.  white 
matter  in  Figure  1).  The  optical  fluence  measured  by  the  probe 
was  in  turn,  calibrated  using  a  calorimeter  (Gold  standard)  and 
showed  a  linear  and  isotropic  response  to  input  NIR  fluence  at 
different  incidence  angles  (<5%  variation).  Our  MC 
simulation  studies  showed  that  99%  of  photons  are  attenuated 
by  the  bone  (1cm),  which  was  confirmed  by  probe 
measurements  in  an  empty  goat  skull.  The  segmentation  code 
(with  iterative  Monte  Carlo)  can  be  used  to  identify  regions  of 
least  bone  thickness  and  determine  the  path  of  least 
attenuation  to  couple  maximum  power  to  brain  tissues.  This 
study  proved  that  an  optimized  Monte  Carlo  could  be 
successfully  used  to  quantify  photon  propagation  in  the  brain 
and  thus  predict  the  rate  of  drug  release  in  target  tissues  to 
depths  up  to  5cm  in  brain  tissues. 
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Fig.  2.  Monte  Carlo  based  treatment  protocol.  Head  CT  image  densities  (in 
Hounsfield  Units/HU)  are  segmented  and  translated  into  optical  properties  of 
the  brain  tissue  (white/gray  matter,  CSF,  skull  bone,  etc.).  Monte  Carlo 
routine  is  then  used  to  simulate  light  propagation  through  the  brain.  [CT 
image  reference:  http://www.mpoullis.com/thorcd/rad/23.htm]16 
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Our  preliminary  longitudinal  study  has  indicated  that  DCE-CT 
is  an  effective  tool  to  study  tumor  progression  and  physiology 
(Fig  3)  of  glioma  tumors.  We  propose  to  monitor  these  effects 
by  studying  the  change  in  fractional  plasma  volume  and 
fractional  interstitial  volume,  both  of  which  are  excellent 
indicators  of  vascular  physiology  (Fig  3). 
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Fig.  3.  Monitoring  cancer  therapy  through  Dynamic  Contrast  Enhanced  CT 
(DCE-CT).  (A)  Schematic  of  two  compartmental  model.12  (B)  Parametric 
images  of  glioblastoma  (orthotopic  model)  grown  in  mice.  In  this  longitudinal 
study,  parametric  images  were  obtained  of  perfusion,  permeability,  fractional 
plasma  volume  and  fractional  interstitial  volume  of  glioma  tumors  grown  in 
mice  using  the  two  compartmental  model.  These  parameters  can  be  used  to 
quantify  the  efficacy  of  treatment,  thus  proving  effectiveness  of  using  DCE- 
CT  to  monitor  targeted  dmg  therapy. 


IV.  Discussion 

In  this  study,  we  have  assumed  that  the  quantity  of  dmg 
release  is  proportional  to  the  optical  power  coupled  to  the 
tissues  (as  shown  in  previous  in-vitro  studies).  Future  studies 
will  focus  on  quantifying  the  optical  dmg  release  threshold  by 
determining  the  relationship  between  optical  power  coupled 
and  amount  of  drug  released.  Further  studies  will  involve 
adjusting  beam  profiles  to  optimize  the  coupling  and 
localization  of  optical  energy  in  the  brain  for  maximal  dmg 
release.  In  order  to  monitor  tumor  response  to  therapy,  DCE- 
CT  will  be  used  to  further  study  vascular  normalization17  in 
tumors  caused  by  anti-angiogenic  effects  of  lapatinib  dmg, 
which  can  help  in  optimizing  therapeutic  efficacy. 


V.  Conclusion 

This  study  has  successfully  demonstrated  the  potential  of 
using  a  fast  voxel  based  Monte  Carlo  based  routine  to  predict 
the  rate  of  dmg  release  in  heterogeneous  brain  tissue  by  using 
Computed  Tomography  images  to  segment  and  translate  tissue 
boundaries  of  the  brain.  An  optimized  therapeutic  protocol 
with  a  fast  Monte  Carlo  routine  can  be  used  to  determine  the 
best  irradiation  conditions  for  targeted  dmg  release  in  brain 
tissue  within  clinically  relevant  time  frames. 
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ABSTRACT 

Purpose:  The  purpose  of  this  study  is  to  develop  an  alternate  empirical  approach  to  estimate  near-infra-red  (NIR)  photon 
propagation  and  quantify  optically  induced  drug  release  in  brain  metastasis,  without  relying  on  computationally 
expensive  Monte  Carlo  techniques  (gold  standard).  Targeted  drug  delivery  with  optically  induced  drug  release  is  a  non- 
invasive  means  to  treat  cancers  and  metastasis.  This  study  is  part  of  a  larger  project  to  treat  brain  metastasis  by 
delivering  lapatinib-drug-nanocomplexes  and  activating  NIR-induced  drug  release.  The  empirical  model  was  developed 
using  a  weighted  approach  to  estimate  photon  scattering  in  tissues  and  calibrated  using  a  GPU  based  3D  Monte  Carlo. 
The  empirical  model  was  developed  and  tested  against  Monte  Carlo  in  optical  brain  phantoms  for  pencil  beams  (width 
1mm)  and  broad  beams  (width  10mm).  The  empirical  algorithm  was  tested  against  the  Monte  Carlo  for  different  albedos 
along  with  diffusion  equation  and  in  simulated  brain  phantoms  resembling  white -matter  (p*’ =8. 25  mm"1,  pa=0. 005mm'1) 
and  gray-matter  (ps’=2.45mm1,  pa=0. 035mm'1)  at  wavelength  800nm.  The  goodness  of  fit  between  the  two  models  was 
determined  using  coefficient  of  determination  (R- squared  analysis).  Preliminary  results  show  the  Empirical  algorithm 
matches  Monte  Carlo  simulated  fluence  over  a  wide  range  of  albedo  (0.7  to  0.99),  while  the  diffusion  equation  fails  for 
lower  albedo.  The  photon  fluence  generated  by  empirical  code  matched  the  Monte  Carlo  in  homogeneous  phantoms 
(R2=0.99).  While  GPU  based  Monte  Carlo  achieved  300X  acceleration  compared  to  earlier  CPU  based  models,  the 
empirical  code  is  700X  faster  than  the  Monte  Carlo  for  a  typical  super-Gaussian  laser  beam. 

Keywords:  Empirical  algorithm  for  light  propagation,  NIR  photon  propagation,  Monte  Carlo  simulation,  Optical 
therapy,  Optically  induced  drug  release,  Tissue  Optics,  Light  Scattering 


1.  INTRODUCTION 

Optical  therapy  using  NIR  light  delivery  requires  the  estimation  of  photon  fluence  in  tissues  and  heterogeneous  media. 
This  study  is  a  part  of  a  larger  project  which  seeks  to  selectively  deliver  macrophages  laden  with  lapatinib -gold- 
nanocomplex  drug  molecules  to  brain  metastasis  with  over-expressed  HER2+  receptors  and  activate  targeted  drug 
release  using  NIR  photons.1,2  The  effectiveness  of  achieving  a  targeted  drug  delivery  and  optically  activated  drug  release 
with  a  Monte  Carlo  based  optical  therapeutic  protocol  has  been  demonstrated  in  previous  studies.1'5  The  design  of  an 
optical  therapeutic  protocol  requires  the  evaluation  of  various  illumination  conditions  such  as  the  angles  of  incidence  and 
source  positions  for  optimum  drug  release  in  the  metastatic  lesions  of  the  brain.5  While  the  Monte  Carlo  is  the  gold 
standard  for  photon  estimation  in  3D  media,  it  is  computationally  expensive  to  be  used  over  a  large  number  of  iterations. 
Hence  there  is  a  need  for  a  computationally  efficient  evaluation  algorithm  which  can  reduce  the  possible  illumination 
conditions  to  a  few  iterations,  while  minimizing  error  in  fluence  estimation.  This  algorithm  should  also  be  able  to 
interface  seamlessly  with  any  imaging  modality  system  to  effectively  localize  fluence  distributions  with  respect  to 
anatomical  and  physiological  information  provided  by  imaging  modalities  such  as  CT,  Ultrasound,  PET,  MRI,  etc.  In  the 
current  work  we  have  developed  a  fast  voxel  based  empirical  approach  which  can  estimate  photon  fluence  in  the  brain 
and  evaluate  different  illuminations  before  using  Monte  Carlo  as  a  final  tool  for  validation.  Being  a  3D  voxel  based 
method;  this  algorithm  can  be  easily  integrated  with  other  imaging  modalities  to  evaluate  fluence  distributions  in  the 
brain  and  other  tissues,  by  co -registering  3D  fluence  images  with  images  from  CT,  MR,  Ultrasound,  etc. 

Light  propagation  studies  in  3D  media  rely  on  the  use  of  Monte  Carlo  simulations  as  the  gold  standard,  to  estimate  the 
photon  fluence  distribution  in  three  dimensional  units  of  volume  known  as  voxels.  The  Monte  Carlo  model  uses 
probability  distribution  to  generate  a  random  walk  of  photons  in  three  dimensional  media  and  relies  on  generating  a  large 
number  of  events  to  produce  statistically  significant  fluence  distributions  in  optical  media  with  specific  absorption  and 
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scattering  properties.6'10  The  photon  fiuence  deposited  (and  traversing)  within  each  voxel  depends  on  the  specified 
optical  properties  assigned  by  the  user,  such  as  the  absorption  coefficient  (pa),  scattering  coefficient  (ps),  anisotropy 
factor  (g)  and  refractive  index  (n).  The  traditional  Monte  Carlo  software  is  computationally  expensive  as  it  relies  on 
generating  millions  of  photons  and  tracking  each  photon  as  it  propagates  through  the  medium.  Recent  research  in  GPU 
based  Monte  Carlo  by  the  Fang  group  has  resulted  in  speeds  up  to  300X  faster  than  the  traditional  CPU  based  versions.6 
However  for  large  optical  datasets,  such  as  those  obtained  and  translated  from  CT  or  MR  images,  and  under  wide  beam 
illumination  conditions,  the  Monte  Carlo  suffers  from  slow  computation  speeds  due  to  large  photon  propagation  times  in 
highly  scattering  tissues.  This  results  in  a  drawback  and  slowdown  in  the  process  of  planning  targeted  light  delivery  to 
specific  regions  within  the  human  body.  In  this  study  we  have  designed  an  alternate  empirical  approach  to  light 
propagation  which  can  estimate  the  photon  fiuence  in  tissues  in  complex  three  dimensional  media  to  provide  a  fast 
preliminary  approach  to  light  estimation  while  planning  targeted  NIR  light  delivery.  The  purpose  of  this  study  is  to 
devise  a  fast  means  to  estimate  NIR  photon  fiuence  under  different  illumination  conditions  in  an  optical  treatment  plan 
before  using  Monte  Carlo  for  final  verification  and  corrections. 

The  empirical  approach  estimates  the  fiuence  distribution  in  a  voxelated  medium  by  using  a  weighted  approach  to 
fiuence  propagation  from  one  layer  to  the  next.  The  assumption  used  is  that  the  photon  fiuence  in  a  voxel  is  a  weighted 
sum  of  the  fluences  in  the  neighboring  voxels  of  the  preceding  layer.  The  weights  are  arbitrarily  assigned  to  each  voxel 
based  on  its  optical  properties,  mainly  the  scattering  coefficient.  The  accuracy  of  the  empirical  approach  depends  on  the 
optimization  method  used  to  estimate  the  weights  for  each  set  of  optical  properties.  The  empirical  approach  simulates 
photon  scatter  by  dividing  the  fiuence  into  scattering  vectors  in  the  forward  and  diagonal  directions,  while  photon 
absorption  is  simulated  using  the  Beer  Lambert’s  rule.11  An  optimization  routine  is  used  to  estimate  the  weights  assigned 
to  each  voxel  by  comparing  with  a  pencil  beam  simulated  by  the  Monte  Carlo.  In  this  study  we  estimated  the  weights  for 
typical  brain  tissues  such  as  Gray  Matter,  White  Matter,  Skull  bone  and  Astrocytoma  with  known  optical  properties.12'14 
These  weights  were  then  tested  in  homogeneous  and  heterogeneous  brain  phantoms  under  pencil  beam  and  broad  beam 
illumination.  The  Monte  Carlo  has  been  used  as  the  gold  standard  for  evaluating  the  performance  of  this  algorithm.  The 
3D  GPU  based  Monte  Carlo,  developed  by  the  Fang  group,  has  an  acceleration  of  300X  over  the  traditional  CPU  based 
version6,  while  the  Empirical  approach  was  seen  to  be  700X  faster  than  the  GPU  based  Monte  Carlo  for  a  super- 
Gaussian  beam. 


2.  MATERIALS  AND  METHODS 


2.1  Algorithm  Design 

The  empirical  approach  approximates  the  photon  scatter  from  one  voxel  to  the  next  by  assigning  weights  to  the 
scattering  directions.  The  main  goal  of  this  algorithm  is  to  approximate  the  total  scattering  events  within  a  voxel  into  14 
component  vectors.  Since  tissues  are  mostly  forward  propagating  with  anisotropies  approximately  equal  to  or  greater 
than  0.9,  it  can  be  assumed  that  most  photon  scatter  directions  can  be  approximated  into  forward  and  diagonal  vectors. 
Hence,  each  scattering  event  is  divided  into  4  scattering  components:  one  forward  component  with  the  wtl  and  4 
diagonal  components  with  weight  wt2  as  shown  in  the  figures.  The  fiuence  enters  a  voxel  from  14  neighboring  voxels. 
Each  incoming  fiuence  component  is  then  scattered  into  one  forward  component  and  4  diagonal  components.  The  total 
fiuence  in  each  individual  direction  is  accumulated  and  absorbed  within  a  voxel  based  on  the  absorption  coefficient.  The 
remaining  fiuence  then  exits  the  voxel  into  the  14  neighboring  voxels.  It  is  important  to  note  that  the  fiuence  entering  a 
voxel  is  equal  to  the  sum  of  fiuence  absorbed  within  the  voxel  and  the  fiuence  exiting  the  voxel.  This  is  achieved  by 
ensuring  that  the  sum  of  the  weights  used  is  equal  to  one  [i.e.  (wtl  +  4  x  wt2)  =1]. 

The  steps  of  the  algorithm  can  be  summarized  as  follows: 

1)  Assign  the  incident  photon  fiuence  (e.g.  Gaussian  beam)  to  the  first  layer  of  voxels.  Here  we  have  assumed  that  the 
photons  enter  the  volume  (of  dimensions  dimlxdim2xdim3)  at  z=l  (first  layer). 

2)  Begin  the  propagation  from  layer  z=l  to  z=dim3. 
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3)  For  every  voxel  in  a  layer  z=k: 

a)  Calculate  the  14  directional  components  from  the  incident  photon  fluence.  For  the  first  layer  the  fluence  is 
assumed  to  be  in  the  forward  direction.  The  fluence  from  neighboring  voxels  is  accumulated  in  each  direction 
as  follows: 

Forward  Component  =  If(x,y,z)  =  £  (forward  components  x  wtl)  +  X(forward  diagonal  components  x  wt2) 

Forward  Diagonal  Component  =  Id(x,y,z)  =  X  (forward  components  x  wt2)  +  X(  forward  diagonal  components 
x  wtl)  +  +  X(side  components  x  wt2) 

Side  Component  =  Is(x,y,z)  =  £(  forward  diagonal  components  x  wt2)  +  £(side  components  x  wtl) 

+  X(  backward  diagonal  components  x  wt2) 

Backward  Diagonal  Component  =  Idb(x,y,z)  =  JXside  components  x  wt2)  +  ^(backward  diagonal  components 
xwtl)  +  I(  backward  component  x  wt2) 

Backward  Component  =  Ib(x,y,z)  =  X  (backward  component  x  wtl)  +  ^(backward  diagonal  components  x  wt2) 
Total  Fluence  =  I(x,y,z)  =  If(x,y,z)  +  XId(x,y,z)  +  JTs(x,y,z)  +  XIdb(x,y,z)  +  Ib(x,y,z) 

The  forward  component  (If)  refer  to  the  direction  of  propagation  of  the  main  fluence  beam  (or  ‘forward 
scatter’),  while  the  backward  component  (lb)  refers  to  the  180  degree  scatter  (or  backscatter)  which  is  opposite 
to  the  incident  photon  fluence  vector.  The  diagonal  components  are  divided  into  4  forward  diagonals  (Id)  and  4 
backward  diagonals  (Idb).  Thus  each  direction  is  associated  with  a  neighboring  voxel.  The  total  fluence  within 
a  voxel  is  obtained  by  summing  up  the  fluences  in  each  component  directions.  Note  that  the  above  directions 
can  be  arbitrarily  chosen  according  to  the  experimental  setup  or  at  the  convenience  of  the  user. 

b)  Simulate  photon  absorption  by  depositing  the  photon  fluence  into  the  voxel  based  on  its  absorption  coefficient 

(Ma): 

Absorbed  Fluence  =  Iabs(x,y,z)  =  I(x,y,z)  x  { l-exp(-pa  x  dl) } 

(where  dl  =  path  length) 

4)  Once  fluences  have  been  accumulated  over  the  current  layer,  increment  z  to  proceed  to  the  next  layer.  Thus  by 
repeating  this  process  for  all  layers  and  updating  the  directional  components  of  fluences  for  every  voxel,  the 
empirical  approach  propagates  the  fluence  from  one  layer  to  the  next  across  the  volume. 

2.2  Optimization  Phase  (Pencil  beam  simulations) 

Optical  properties  of  the  medium  (scattering  coefficient  ps,  anisotropy  factor  g,  absorption  coefficient  pa,  albedo  a’) 
define  the  photon  fluence  generated  by  the  empirical  approach.  The  weights  assigned  to  each  voxel  depend  on  its 
scattering  property,  while  the  absorption  coefficient  determines  the  amount  of  fluence  deposited  and  energy  absorbed. 
The  weights  were  derived  during  the  optimization  phase  by  comparing  the  exponential  attenuation  profiles  of  the 
Empirical  approach  to  the  Monte  Carlo  using  root  mean  square  error  and  coefficient  of  correlation  criteria.  Previous 
studies  have  shown  the  applications  as  well  as  the  limitations  of  the  diffusion  equation  in  estimating  photon  fluence  for 
media  with  low  albedos.10,15'22  In  this  study  we  first  tested  the  performance  of  the  Empirical  approach  with  respect  to  the 
Monte  Carlo  and  the  diffusion  equation  for  media  with  different  albedos  from  0.7  to  0.99.  The  albedo  (a’)10  is  given  by 
the  equation: 


(i-gK 

(l-<?Vs  +  Ma 


(1) 


Next,  the  empirical  algorithm  was  trained  to  derive  the  forward  (wtl)  and  diagonal  (wt2)  weights  for  brain  tissues  such 
as  white  matter  (p^O.OOSmm1,  ps=55rnrrf g=0.85,  a’=0.9994),  gray  matter  (p^O.OSSmm1,  ps=70mm'',  g=0.965, 
a’=0.9859),  skull  bone  (pa=0. 024mm'1,  p^lS^mm'1,  g=0.9,  a’=0.9871)  and  tumor  (astrocytoma)  (p^O.Smm1, 
ps=500mm_1,  g=0.88,  a’=0.995).12'14  A  pencil  beam  was  simulated  using  the  Monte  Carlo  using  the  known  optical 


25 


properties  of  each  of  the  above  brain  tissues.12'14  The  Empirical  algorithm  was  then  trained  to  emulate  the  same 
exponential  decay  profile  by  incrementing  the  forward  scattering  weight  from  0  to  1,  while  the  diagonal  was  adjusted  to 
ensure  that  wtl  +  (4  x  wt2)  =1.  The  weights  derived  from  the  optimization  phase  were  used  for  the  broad  beam 
simulations  in  homogeneous  and  heterogeneous  phantoms. 

2.3  Broad  beam  simulations  in  Homogeneous  and  Heterogeneous  Brain  phantoms 

The  empirical  weights  obtained  from  the  pencil  beam  training  phase  were  used  for  further  simulations  using  broad 
beams  in  homogeneous  and  heterogeneous  brain  phantoms.  The  Gaussian  beam  generated  by  the  Monte  Carlo  was 
matched  to  that  of  the  Empirical  algorithm.  An  optical  property  map  was  generated  in  which  every  voxel  was  assigned  a 
specific  weight  and  absorption  coefficient  based  on  the  tissue  property.  The  Gaussian  beam  distribution  was  first  tested 
in  homogenous  phantoms  of  each  of  the  brain  tissues:  skull  bone,  gray  matter,  white  matter  and  tumor  (astrocytoma). 
The  weights  derived  from  the  pencil  beam  optimization  phase  were  used  in  the  broad  beam  simulation  to  test  their 
validity.  Next,  a  software  phantom  of  a  human  head  with  15  cm  diameter  was  created  using  concentric  spheres  to  denote 
the  different  tissue  types.  The  incident  beam  had  a  Gaussian  profile.  The  fluence  profiles  and  absorption  energy  maps 
generated  Monte  Carlo  and  Empirical  approach  were  compared.  The  execution  times  for  generation  of  a  broad  beam 
profile  was  compared  between  the  two  models. 


3.  RESULTS 


3.1  Pencil  beam  simulations  (Optimization  Phase) 


In  this  study,  the  empirical  software  was  trained  in  an  iterative  process  by  comparing  its  exponential  decay  profiles  to 
that  of  the  Monte  Carlo  using  pencil  beams  in  optically  homogeneous  media.  The  Empirical  approach  had  a  better  fit  to 
the  Monte  Carlo  as  compared  to  the  diffusion  equation,  especially  for  low  albedos.  The  weights  wtl  and  wt2  were 
derived  and  optimized  for  each  of  the  brain  tissues  (white  matter,  gray  matter,  skull  bone  and  tumor)  by  using  the  profile 
which  best  matched  the  Monte  Carlo.  The  absorption  coefficient  assigned  was  the  same  as  has  been  obtained  from 
previous  studies.  Figure  3  shows  the  pencil  beam  simulation  profiles  obtained  by  the  Monte  Carlo  and  the  Empirical 
algorithm. 
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Figure  3.  (a)  Training  Phase.  The  figure  shows  pencil  beam  profiles  along  z  axis  between  the  Monte  Carlo  and  the 
Empirical  algorithm.  The  empirical  algorithm  is  trained  to  emulate  the  same  exponential  fluence  attenuation  as  that  of 
the  Monte  Carlo  by  adjusting  and  optimizing  the  weights  wtl  and  wt2.  (b)  Broad  beam  simulations:  The  empirical 
algorithm  was  compared  to  the  Monte  Carlo  for  incident  Gaussian  beam  profiles.  Overall  the  fluence  profile  of  the 
Empirical  closely  followed  that  of  the  Monte  Carlo. 
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Figure  3(b)  compares  the  Gaussian  beam  attenuation  along  beam  axis  between  the  Empirical  and  the  Monte  Carlo.  The 
Gaussian  beam  analysis  shows  that  the  empirical  approach  can  emulate  the  Monte  Carlo  fluence  for  broad  beam 
simulations.  The  acceleration  provided  by  the  Empirical  algorithm  was  approximately  700X  as  compared  to  the  Monte 
Carlo  for  a  broad  beam.  This  however  depends  on  the  tissue  optical  properties  as  the  Monte  Carlo  simulations  take 
longer  execution  times  for  high  scattering  tissues  and  are  relatively  fast  for  low  scattering  media.  The  differences 
between  the  Monte  Carlo  and  the  Empirical  approaches  can  be  significant  at  larger  depths  and  hence  need  to  be 
accounted  for. 
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3.3  Heterogeneous  Brain  phantom  simulations  (Gaussian  beam) 

Figure  5  shows  the  fluence  profiles  and  absorption  energy  maps  generated  Monte  Carlo  and  Empirical  approach  for  the 
head  phantom.  The  errors/differences  induced  by  the  tissue  heterogeneity  are  seen  between  the  two  models  and  needs  to 
be  studied  further.  However  it  is  seen  that  the  exponential  profile  of  the  Empirical  more  or  less  follows  that  of  the  Monte 
Carlo.  In  the  absorbed  energy  images  generated  in  Figure  5  (b)  and  (c),  the  profile  of  the  tumor  absorption  is  clearly 
visible  in  the  Empirical  and  can  be  demarcated  within  the  Monte  Carlo.  The  Monte  Carlo  also  generates  a  large  amount 
of  statistical  noise  at  larger  depths  which  is  not  seen  in  the  Empirical  model.  The  time  for  execution  of  the  Monte  Carlo 
for  a  broad  beam  in  the  heterogeneous  head  phantom  is  3.01  seconds  as  compared  to  9  hours  and  47  minutes  required  by 
the  Monte  Carlo.  Thus  the  Empirical  model  gives  a  smoother  approximation  of  the  fluence  profile  as  compared  to  the 
Monte  Carlo  at  a  fraction  of  the  time  (1170 IX  acceleration)  and  thus  could  provide  an  alternative  to  the  Monte  Carlo  in 
cases  where  an  approximate  idea  of  fluence  distribution  is  needed  while  planning  optical  therapy. 


4.  DISCUSSION 


4.1  Implications  and  Future  Research 

The  current  study  has  looked  at  providing  an  alternate  voxel  based  approach  of  estimating  photon  fluence  in  3D 
heterogeneous  tissues  besides  the  diffusion  approximation  and  Monte  Carlo  methods.  The  Monte  Carlo  is  by  far  the 
most  accurate  method  for  simulating  photon  deposition  in  tissues.  However  due  to  the  large  execution  times,  even  after 
being  implemented  on  parallel  processing  platforms,  there  is  a  need  for  a  method  for  fast  preliminary  estimation  of 
fluence  within  a  certain  error  limit.  The  empirical  approach  is  one  such  method  which  is  fast  and  can  approximate  the 
fluence  in  heterogeneous  tissues.  By  generating  fluence  in  a  3D  voxelated  format,  it  can  be  easily  integrated  with 
multiple  imaging  modalities.  More  study  is  however  necessary  to  understand  the  errors  induced  by  the  Empirical 
algorithm  at  large  depths  as  seen  in  the  pencil  beam  and  broad  beam  fluence  profiles.  These  errors  seen  can  be  large  at 
depths  greater  than  2cm  in  heterogeneous  phantoms.  Hence  the  empirical  approach,  at  this  time,  cannot  be  used  for 
absolute  quantification  of  photon  fluence.  It  can  however  give  a  good  approximation  of  the  relative  fluence  distribution 
in  tissues  with  respect  to  the  input  photon  fluence.  The  errors  in  quantification  can  be  reduced  by  further  research  into 
optimization  methods  to  match  the  Monte  Carlo.  The  training  methodology  used  to  derive  the  weights  could  be 
improved  to  account  for  the  differences  in  the  fluence  profiles.  Iterative  routines  and  algorithms  with  merit  functions  can 
significantly  enhance  the  accuracy  of  the  Empirical  approach. 

4.2  Conclusion 

The  ability  of  the  Empirical  approach  to  approximately  emulate  the  fluence  distribution  of  the  Monte  Carlo  has  been 
studied  in  homogeneous  and  heterogeneous  media.  Even  though  the  empirical  algorithm  cannot  be  used  for  absolute 
quantification,  it  can  be  used  in  optical  therapy  to  evaluate  the  merits  of  different  angles  of  illuminations  and  source 
positions  which  is  a  prerequisite  for  determining  the  placement  of  an  optical  source.  The  3D  GPU  based  Monte  Carlo 
can  then  be  used  for  absolute  quantification  and  validation  of  photon  energy  distribution.  The  acceleration  provided  by 
the  Empirical  algorithm  can  be  used  to  reduce  the  iterations  of  the  Monte  Carlo  by  identifying  only  a  few  possible 
illumination  possibilities.  Thus  the  Empirical  algorithm  has  the  potential  to  be  used  as  a  fast  preliminary  tool  to  devise 
optical  therapies  by  estimating  relative  fluence  distribution  in  tissues  prior  to  validation  and  correction  by  using  Monte 
Carlo. 
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A  three  dimensional  quantification  of  photon  energy  distribution  in  tissues  is  critical  in  designing  optical 
therapeutic  protocols  to  trigger  light  activated  drug  release  in  tissues.  The  objective  of  this  study  is  to 
develop  software  and  dosimetry  tools  to  simulate,  characterize  and  cross-validate  the  local  photon  fluence 
in  heterogeneous  brain  tissue,  with  an  aim  to  quantify  the  effects  of  photo-activated  drug  release  in  brain 
tumors.  A  GPU  enhanced  voxel-based  Monte  Carlo  code  was  developed  to  simulate  the  transport  of  near 
infrared  photons  with  differing  laser  beam  profiles  (i.e. ,  flat  and  super-gaussian)  within  phantoms  of  white 
and  gray  matter.  From  these  simulations,  the  local  photon  fluence  and  tissue  dosimetric  distribution  was 
validated  through  the  implementation  of  an  optical  dosimetry  probe  with  an  isotropic  acceptance  and  1.5mm 
diameter.  To  translate  these  concepts  into  clinical  practice,  empirical  models  of  photon  transport  were 
devised,  calibrated  to  MC,  and  realized  in  CUDA  to  provide  3D  fluence  and  optical  dosimetric  maps  in 
real-time.  These  tools  provide  the  capability  to  develop  and  optimize  treatment  plans  for  optimal  release 
of  pharmaceuticals  to  metastatic  breast  cancer  in  the  brain.  Future  work  will  test  and  validate  these  novel 
delivery  and  release  mechanisms  in  vivo. 


1.  Introduction  Monte  Carlo  software,  enhanced  with  a  quick  prelimi- 

Optical  therapy  has  been  used  extensively  to  treat  can-  nary  empirical  photon  propagation  model  can  be  used 
cer,  e.g.  directly  in  Photo  Dynamic  Therapy  (PDT)  to  devise  an  optimal  therapeutic  plan  within  a  clinically 

[1]  and  indirectly  to  trigger  release  of  drug  molecules  relevant  time  frame.  In  this  study  we  have  investigated 

from  nanocomplexes  in  optically  targeted  therapy  [2,  3].  this  approach  and  validated  the  Monte  Carlo  and  the 
Near-infra-red  (NIR)  light  source  is  preferred  for  deep  preliminary  empirical  model  using  a  specially  developed 

tissue  activation  due  to  its  high  penetration  and  low  ’’linear  and  isotropic”  optical  dosimetry  probe  in  brain 

optical  attenuation  in  tissues.  This  study  is  part  of  a  tissue  phantoms.  The  optical  dosimetry  probe  served  as 

project  which  seeks  to  develop  optically  activated  drug  the  gold  standard  for  validating  the  Monte  Carlo  model 

release  in  the  brain  and  has  successfully  demonstrated  and  has  a  better  design  compared  to  previous  models, 

a  macrophage  based  Trojan  Horse  delivery  and  activa¬ 
tion  of  molecules  into  brain  metastasis  [2,  3].  The  broad  Monte  Carlo  simulation  has  been  used  as  the  gold 

goal  is  to  develop  an  optically  targeted  therapy  for  brain  standard  for  photon  propagation  studies  in  heteroge- 

metastasis  by  delivering  lapatinib-gold  nanocomplex-  neous  media  and  uses  probability  distribution  to  gen- 

laden-macrophages  to  metastastic  lesions  and  optically  erate  the  random  walk  of  photons  based  on  the  absorp- 

triggering  drug  release,  thus  combining  the  effects  of  tar-  tion  and  scattering  properties  of  tissues  [4-8] .  Tradi- 

geted  drug  delivery  with  targeted  activation.  With  the  tional  CPU  based  Monte  Carlo  models  suffer  from  slow 

rate/quantity  of  drug  release  being  directly  dependent  computation  speeds,  which  makes  them  non- viable  for 

on  the  optical  power  delivered  [2,  3],  a  3D  Monte  Carlo  application  in  clinical  settings,  especially  in  high  scat- 

(MC)  based  protocol  can  be  used  to  estimate/optimize  tering  tissues  (e.g.  white  matter).  The  GPU  based 

the  photon  energy  deposition  (and  subsequent  drug  ac-  3D  Monte  Carlo  (MCX  model)  used  in  this  study  was 

tivation)  in  brain  tissues  [10].  A  fast  GPU  based  3D  developed  by  the  Fang  group  [9]  and  has  been  shown 

_  to  have  a  300X  acceleration  compared  to  the  previous 

CPU  based  models.  The  CUDA  acceleration  provided 

*  Corresponding  author:  aprabhuv@purdue.edu  by  this  model  helps  to  achieve  complete  brain  simulation 

t  Corresponding  author:  fangq@nmr.mgh.harvard.edu  of  a  broad  beam  NIR  source  within  minutes,  as  COm- 

*  Corresponding  author:  susan.clare@northwestern.edu  pared  to  hours  of  simulation  time  using  CPU  based  mod- 

§  Corresponding  author:  kstantz@purdue.edu  els  jgj  Furthermore  the  model  provides  a  voxel  based 
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3D  distribution  of  photon  fluence  and  energy  deposited, 
which  can  be  co-registered  with  different  medical  imag¬ 
ing  modalities  (e.g.  CT,  MR,  Ultrasound,  etc),  thus  im¬ 
proving  its  clinical  viability  for  diagnosis  and  therapy. 
In  order  to  assess  the  effects  of  different  illumination 
conditions  on  NIR  photon  distribution  (and  subsequent 
drug  release)  in  brain  tissues,  we  customized  the  Monte 
Carlo  to  emulate  these  user  defined  broad  beam  profiles: 
(a)  diffuse  broad  beam  emanating  from  an  integrating 
sphere  and  (b)  super- Gaussian  beam  distribution  from 
an  optical  fiberguide  connected  to  a  pulsating  NIR  laser 
(wavelength  800nm).  Optical  brain  phantoms  resem¬ 
bling  white  and  gray  matter  were  designed  using  In¬ 
dia  ink  (absorber)  and  Intralipid  (scatterer)  to  resemble 
the  absorption  and  scattering  properties  of  brain  tissues. 
The  GPU  based  Monte  Carlo  was  validated  by  compar¬ 
ing  the  MC  generated  fluence  to  localized  fluence  distri¬ 
bution  measured  using  optical  dosimetry  probes.  The 
use  of  a  fast  3D  Monte  Carlo  allows  the  development  of 
an  optical  therapeutic  protocol  to  predict  photon  dis¬ 
tribution  and  targeted  drug  release  in  the  brain  in  a 
clinically  relevant  timeframe  [10]. 

Optical  dosimetry  is  a  technique  of  measuring  local¬ 
ized  photon  distribution  in  tissues.  Optical  dosimetry 
probes  with  an  isotropic  response  have  been  extensively 
used  in  photodynamic  therapy  (PDT)  to  measure  pho¬ 
ton  fluence  as  well  as  to  deliver  light  locally  to  tissues 
[1,  11,  13-16].  The  small  size  of  the  dosimetry  probe 
(1.5mm  diameter)  allows  it  to  be  used  in  localized  light 
measurement  in  tissues  as  well  as  optical  phantoms,  un¬ 
like  standard  optical  fluence  measurement  devices  such 
as  calorimeters,  which  cannot  measure  isotropic  fluence 
and  are  limited  by  their  size.  In  this  study,  we  designed 
two  types  of  optical  dosimetry  probes,  with  tips  made 
of  (a)  titanium  and  (b)  nylon,  and  evaluated  their  ef¬ 
fectiveness  (probe  response,  linearity  and  isotropicity) 
in  measuring  localized  fluence  in  tissue  phantoms.  The 
titanium  based  dosimetry  probes  were  found  to  be  su¬ 
perior  in  performance,  with  a  linear  and  isotropic  re¬ 
sponse,  compared  to  traditional  Nylon  based  probes. 
The  optical  dosimetry  probes  were  used  to  quantify  the 
near  infrared  (NIR)  photon  propagation  in  brain  tissue 
phantoms  and  validate  the  3D  Monte  Carlo  model.  The 
dosimetry  probes  thus  helped  in  calibration  of  the  Monte 
Carlo  software  as  a  simulation  tool  to  accurately  predict 
NIR  photon  energy  deposition  in  brain  tissues  and  to 
emulate  various  light  sources  for  clinical  applications  in 
optical  therapy. 

In  addition  to  the  Monte  Carlo  model  used  in  this 
study,  we  designed  an  empirical  photon  propagation 
model  to  estimate  photon  fluence  in  tissues  and  het¬ 
erogeneous  media  at  a  fraction  of  time  taken  by  the 
GPU  based  Monte  Carlo  [17].  This  Empirical  approach 
is  based  on  the  photon  convolution  method  and  is 
a  weight  based  photon  propagation  algorithm.  This 
model  assumes  that  the  photon  fluence  in  a  voxel  is  the 
weighted  sum  of  the  fluences  in  the  neighboring  voxels 
in  the  preceding  layer.  This  approach  approximates 


photon  scatter  between  one  layer  to  the  next  using 
arbitrarily  assigned  weights,  and  simulates  absorption 
by  depositing  photon  fluence  in  each  voxel  before 
proceeding  to  the  next  layer  [17].  The  absorption  and 
scattering  weights  used  in  this  method  are  calculated  by 
comparing  it  to  the  Monte  Carlo  model  in  an  ’optimiza¬ 
tion’  phase  (using  R-square  analysis).  We  optimized  and 
tested  the  Empirical  algorithm  against  the  Monte  Carlo 
(gold  standard)  for  typical  brain  tissues  such  as  Gray 
Matter,  White  Matter,  Skull  bone  and  Astrocytoma 
with  known  optical  properties.  The  empirical  approach 
can  be  used  as  a  preliminary  planning  tool  to  narrow 
down  and  optimize  illumination  conditions  for  broad 
beams  in  heterogeneous  brain  tissues  before  using  the 
more  accurate  Monte  Carlo  to  estimate  and  optimize 
photon  energy  distribution  and  subsequent  drug  release. 
A  targeted  MC  based  Optical  Therapeutic  protocol 
enhanced  by  the  Empirical  approach  can  accurately 
predict  and  optimize  photon  activated  drug  release, 
while  sparing  healthy  tissues  in  the  brain  [10]. 

2.  Materials  and  Methods 
2. A.  Optical  Brain  Phantoms 

Optical  brain  phantoms  were  designed  to  emulate  the 
optical  properties  of  the  brain  (absorption  and  reduced 
scattering  coefficients),  namely  white  matter,  gray  mat¬ 
ter  and  skull  bone  at  an  NIR  wavelength  of  800nm 
[18,  19].  Previous  studies  have  shown  that  reliable  op¬ 
tical  phantoms  can  be  built  using  India  ink  (absorber) 
and  intralipid  (scatterer)  [20-22].  India  ink  has  a  negli¬ 
gible  scattering  coefficient  compared  to  intralipid,  while 
absorption  coefficient  of  intralipid  is  in  turn  negligible 
compared  to  India  ink  [20-22].  Thus  the  net  scattering 
coefficient  in  the  phantom  is  contributed  by  intralipid 
while  the  absorption  coefficient  is  due  to  the  absorption 
properties  of  India  ink,  in  addition  to  the  optical  absorp¬ 
tion  by  water  molecules  (pa  =  0.02cm-1)  [23]. 

The  absorption  coefficient  of  India  ink  (Higgins  non 
water  proof  black  ink),  at  800nm  wavelength,  was  deter¬ 
mined  and  validated  using  a  spectro-photometer  (Shi- 
madzu  UVmini  1240  Spectrophotometer),  which  served 
as  the  gold  standard  to  measure  absorption.  In  order  to 
measure  the  scattering  coefficient  of  20%  intralipid  (In¬ 
tralipid  20%  from  Fresenius  Kabi),  a  plexiglass  cuvette 
filled  with  20%  intralipid  (in  various  concentrations)  was 
illuminated  with  a  uniform  NIR  broad  beam  (wavelength 
800nm)  from  an  integrating  sphere,  and  the  attenua¬ 
tion  coefficient  measured  by  translating  the  probe  along 
the  beam  axis  in  the  solution.  The  attenuation  coef¬ 
ficient  measured  by  the  probe  is  due  to  scattering  by 
intralipid  and  absorption  by  water  molecules.  We  then 
matched  the  probe-measured  attenuation  coefficient  us¬ 
ing  Monte  Carlo  (for  similar  illumination  conditions), 
by  adjusting  the  scattering  coefficient  of  the  simulated 
medium.  The  absorption  coefficient  of  the  medium  (by 
water  molecules)  was  kept  constant  at  pa  =  0.02cm-1 
[23].  The  MC  derived  scattering  coefficient  which  gave 
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the  best  match  (highest  R-square  value)  for  the  atten¬ 
uation  coefficient  was  selected  as  the  derived  scattering 
coefficient  value  for  that  concentration. 

Figure  1  shows  the  calibration  curves  of  concentration 
versus  absorbance  of  India  ink  obtained  using  the  photo¬ 
spectrometer,  while  figure  2  shows  the  reduced  scatter¬ 
ing  coefficient  versus  concentration  of  intralipid  obtained 
from  the  Monte  Carlo  method.  These  graphs  show  a  lin¬ 
ear  relationship  between  concentration  and  attenuation 
coefficient  (absorption  in  figure  1  and  scattering  in  figure 
2)  and  was  used  to  design  optical  phantoms  with  specific 
absorption  and  scattering  properties  to  resemble  white 
matter,  gray  matter  and  skull  bone.  It  is  important  to 
note  that  the  linear  relationship  of  concentration  versus 
attenuation  can  change  due  to  miniscule  (yet  important) 
scattering  properties  of  India  ink  and  absorption  prop¬ 
erties  of  intralipid.  These  effects  while  minimal  over  a 
small  range  of  concentrations,  can  be  a  source  of  uncer¬ 
tainties  over  a  larger  range  of  concentrations  in  phantom 
design.  Two  types  of  phantoms  were  built  to  resemble 
white  and  gray  matter.  Previous  studies  by  Van  der  Zee 
[18]  have  determined  the  optical  properties  of  white  and 
gray  matter,  while  the  study  by  Firbank  [19]  determined 
those  of  skull  bone.  These  optical  properties  (absorp¬ 
tion  and  scattering  coefficients  and  anisotropy  factor) 
are  shown  in  table  1  and  were  used  to  design  the  optical 
brain  phantoms  using  the  equations  derived  from  figures 
2  and  3. 


Fig.  1.  Calibration  curve  of  Spectrophotometer  absorbance 
versus  Concentration  of  India  ink. 

In  order  to  determine  the  concentration  of  intralipid 
and  India  ink,  the  linear  relationship  of  concentration 
versus  total  attenuation  was  used. 

2.B.  Dosimetry  probe  design  and  Calibration 

The  optical  dosimetry  probe  consists  of  a  spherical  bulb 
of  highly  scattering  material  attached  to  an  optical  fiber 
[1,  11,  12,  14].  The  optical  fiber  is  connected  to  a  PIN 


Intralipid  Calibration:  Reduced  Scattering  Coefficient  v/s 
Concentration 


Concentration  %  (Intralioid) 


Fig.  2.  Calibration  curve  of  Reduced  scattering  coefficient 
versus  Concentration  of  Intralipid. 


Table  1.  Optical  Properties  of  Brain  Tissues  and  Equivalent 
Phantom  Composition  [18,  19] 


Tissue 

Type 

Absorption 

Coefficient 

pa  (cm-1) 

Scattering 

Coefficient 

ps  (cm-1) 

Anisotropy 

Factor 

g 

Reduced 

Scattering 

Coefficient 

Ps1  (cm-1) 

Intralipid% 

rpJ+2.17481 

Ink% 

r  Pa  —  0 •  0049  -i 

L  11.297  J 

L  15.066  J 

White 

Matter 

0.05 

550 

0.85 

82.5 

7.5 

0 

Gray 

Matter 

0.35 

700 

0.965 

24.5 

2.36 

0.023 

Skull 

Bone 

0.24 

184 

0.9 

18.4 

1.82 

0.0156 

photodiode  circuit  whose  output  voltage  is  read  using  an 
oscilloscope.  The  output  voltage  is  proportional  to  the 
optical  fluence  entering  the  spherical  bulb.  Two  models 
optical  dosimetry  probes,  with  tips  made  of  Nylon  and 
Titanium  dioxide,  were  designed  [26].  The  performance 
of  the  traditionally  used  Nylon  based  probes  [14]  was 
compared  to  that  of  the  novel  Titanium  based  probes  by 
measuring  the  isotropicity  of  response  (in  equatorial  and 
azimuthal  directions)  to  different  angles  of  incidence. 

The  Nylon  spheres  have  been  widely  used  for  photody¬ 
namic  therapy  and  were  made  of  Nylon™  due  to  their 
isotropic  response  [14].  These  were  used  as  the  standard 
of  reference  to  compare  the  performance  of  the  novel  Ti¬ 
tanium  based  probes.  The  titanum  based  probes  have 
the  scattering  tips  made  of  a  mixture  titanium  dioxide 
(Du  Pont  Ti-Pure  R-900)  and  a  clear  two-part  epoxy 
(Tra-Con  BA-F114)  mixed  in  a  ratio  of  9.1  mg  Ti02  to 
1  ml  of  the  epoxy  (mixed)  [26,  27].  The  titanium  probes 
were  selected  since  the  particle  size  of  titanium  dioxide 
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Optical  dosimetry  probe 


Fig.  3.  Setup  for  calibration  of  the  Optical  dosimetry  probe: 
measuring  isotropicity  (in  equatorial  and  azimuthal  planes) 
of  the  probe. 


is  around  410nm  [27],  which  results  in  optimal  scatter¬ 
ing  occurs  820nm,  i.e.  when  the  wavelength  is  approxi¬ 
mately  twice  the  particle  size  [27].  Since  the  wavelength 
used  in  our  studies  is  800nm,  the  titanium  probes  have 
been  designed  to  show  a  higher  sensitivity  than  the  Ny¬ 
lon  probes.  The  titanium  based  probe  tips  were  molded 
within  a  specially  designed  mould  [26]  and  a  fiber  optic 
cable  (BFL48-400,  4mm  diameter,  0.48  numerical  aper¬ 
ture)  was  inserted  within  the  cured  spheres  at  a  depth 
of  approximately  1/3  of  the  sphere’s  radius  [16].  The 
Nylon  based  probes  were  apprximately  3.175mm  diame¬ 
ter  and  were  machined  in  a  multistage  process  [26] .  The 
other  end  of  the  fiber  optic  cable  was  connected  to  PIN 
photodiode  (ET-2030,  Electro  Optics  Technology  Inc.) 
though  a  SMA-905  connnector  and  the  voltage  output 
was  read  using  an  oscilloscope  (Tektronix  TDS3052B). 

Figures  3  shows  the  setup  used  to  test  the  isotrop¬ 
icity  of  the  dosimetry  probes  in  the  equatorial  and  az¬ 
imuthal  directions.  The  probe  tip  is  illuminated  by  laser 
light  emanating  from  an  integrating  sphere  (Melles  Griot 
two  port  integrating  sphere),  which  produces  a  uniform 
broad  beam  of  diffuse  light,  which  is  then  used  to  char¬ 
acterize  the  probe  response.  The  probe  is  rotated  at 
specific  angles  in  the  equatorial  and  azimuthal  planes 
and  its  response  (output  rms  voltage)  is  measured.  The 
center  of  the  probe  tip  was  kept  at  the  same  location 
during  this  rotation,  thus  only  changing  its  angular  ori¬ 
entation.  The  isotropicity  of  the  probe  is  measured  by 
calculating  its  coefficient  of  variation,  which  is  the  ra¬ 
tio  of  standard  deviation  to  mean  of  rms  value  over  all 
the  angles.  The  range  of  angles  were  0  to  360  degrees 
for  equatorial  measurements  and  0  to  150  degrees  for  az¬ 
imuthal  measurements.  The  angular  span  for  Azimuthal 
measurements  was  limited  due  to  physical  limitations  of 


the  setup  (figure  4)  which  prevented  larger  angles  due  to 
presence  of  the  integrating  sphere.  The  individual  rms 
voltage  measured  by  the  probe  was  normalized  to  the 
power  measured  using  a  calorimeter. 

In  the  test  of  linearity  (figure  3(A)),  the  dosimetry 
probe  response  is  validated  using  a  calorimeter  (Ophir 
CE,  [Nova  2  Ophir:  3A-P-SH-VI],  aperture  12mm)  by 
illuminating  with  a  pulsed  laser  beam  of  wavelength  800 
nm,  frequency  10  Hz,  beam  intensity  10.1911  mW/cm2 
and  energy  0.02  J/pulse.  The  probe  and  the  calorime¬ 
ter  are  translated  along  the  beam  axis  (in  air)  and 
the  photon  fluence  is  measured  at  different  positions. 
The  linearity  of  response  (probe  v/s  calorimeter)  of  the 
1.5mm  Titanium  probe  (selected  due  to  its  isotropic¬ 
ity)  is  shown  in  figure  5(B).  The  linearity  of  the  probe 
measurement  (c i)probe ,  mV/mm2)  with  respect  to  the 
calorimeter  measurement  (4>  calorimeter,  mV/mm2)  was 
derived  as  follows  by  fitting  a  line  on  the  measured  data: 

^  _  0 calorimeter  T  1*7296 

Qprobe  ~  590.46  ^  j 

The  probe  voltage  signal  was  filtered  using  a  band¬ 
pass  Butterworth  filter  to  remove  high  frequency  noise 
in  the  probe  response  (figure  6).  The  photon  fluence 
is  proportional  to  integral  of  the  probe  voltage  output 
over  the  time  of  the  laser  pulse  and  was  calculated  by 
summing  over  the  measured  voltages  (in  the  waveform) 
over  the  pulse  duration.  This  is  a  novel  technique  used 
in  this  study  and  was  found  to  be  a  much  better  repre¬ 
sentation  of  the  probe  output  as  compared  to  measuring 
the  absolute  voltage  or  the  root  mean  squared  value. 
The  probe  voltage  (integrated  over  the  pulse  duration) 
and  calorimeter  power  were  normalized  to  the  respective 
surface  areas  [probe  radius  =  1.5mm,  probe  surface  area 
(sphere)  =  7.065mm2,  calorimeter  aperture  =  12mm, 
calorimeter  surface  area  (circle)  =  113.04mm2]. 

2.C.  Absolute  Photon  Quantification  in  tissue-like 
media 

The  linear  relationship  in  equation  1,  used  to  quantify 
the  photon  dose  in  air,  was  modified  to  derive  the  opti¬ 
cal  fluence  in  liquid  brain  phantoms.  The  ratio  of  light 
entering  the  fiber  within  the  probe  depends  on  the  ratio 
of  the  refractive  index  of  the  probe  with  respect  to  the 
medium.  The  refractive  index  of  the  titanium  probe  tip 
is  approximately  2.73  [27]  as  it  is  made  of  rutile  titanium 
dioxide.  Since  the  refractive  index  of  air  is  smaller  than 
that  of  water,  the  probe  response  in  water  is  smaller  as 
more  light  leave  the  probe  in  water  compared  to  air.  The 
reflectance  factor  is  0.41  for  the  probe  to  air  refractive 
index  ratio  of  2.5,  and  0.74  for  the  probe  to  water  ra¬ 
tio  of  1.88  [25].  The  difference  in  the  refractive  index 
and  reflectance  between  probe-air  and  probe-water  in¬ 
terfaces  (water  and  tissues  have  similar  refractive  index 
of  1.329)  must  be  accounted  for  to  quantify  photon  dose 
in  water.  In  order  to  determine  the  calibration  factor  we 
measured  the  probe  response  in  air  and  water  within  a 
cuvette  illuminated  by  a  diffuse  broad  uniform  beam  em- 
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Optical  Dosimetry  Probe 


Plexiglas  Cuvette  filled  with 


NIR  Photon  source: 

1.  Broad  beam  (Integrating 
Sphere) 

2.  Super-Gaussian  beam 
(Optical  Fiberguide) 


>=> 


Fig.  4.  (A) Setup  for  calibration  of  the  Optical  dosimetry 

probe:  measuring  the  linearity  of  the  probe  to  incident  laser 
fluence.  (B) Setup  for  Optical  dosimetry  in  a  liquid  phantom. 
The  optical  dosimetry  probe  is  translated  within  a  plexiglas 
cuvette  filled  with  liquid  phantom,  and  illuminated  by  an 
NIR  source.  This  setup  is  used  for  intralipid  calibration  as 
well  as  for  actual  phantom  studies 


anating  from  an  integrating  sphere  (same  setup  as  figure 
3).  The  probe  was  placed  inside  the  cuvette,  0.2cm  from 
the  surface.  First  the  probe  response  voltage  was  mea¬ 
sured  in  air  (c t>air )  and  then  measured  again  by  filling 
the  cuvette  with  water  ( 4>water ),  without  changing  the 
position  of  the  probe  with  respect  to  the  cuvette  and  the 
laser  source.  By  accounting  for  the  attenuation  by  wa¬ 
ter  (Hwater  —  0.02 cm-1)  for  a  pathlength  (1)  of  0.2cm, 
and  the  beam  spread  (due  to  change  in  refractive  in¬ 
dex)  differences  between  the  two  setups  (cuvette-to-air 
v/s  cuvette-to- water),  the  calibration  factor  (CF)  was 
determined  experimentally  to  be  1.3182,  using  the  fol¬ 
lowing  equation: 


CF  = 


< fiwater 


(2) 


We  neglect  the  attenuation  in  air  in  the  above  cal¬ 
culation,  and  hence  the  probe  response  in  air  at  0.2cm 
from  the  cuvette  surface  is  assumed  to  be  same  as  that 


Fig.  5.  (A)Probe  Isotropicity  of  Titanum  v/s  Nylon  probes. 
The  coefficient  of  variation  (standard  deviation/ average  rms 
voltage)  was  calculated  for  probes  made  of  Nylon  and  Tita¬ 
nium.  The  error  bars  show  indicate  the  minimum  and  max¬ 
imum  values  of  variation  in  individual  probe  response.  The 
Titanium  based  probes  showed  a  more  isotropic  response  in 
the  equatorial  and  azimuthal  directions.  A  1.5mm  diameter 
Titanium  probe  was  selected  for  further  dosimetric  studies. 
(B) Probe  linearity:  graph  of  total  voltage  (area  under  pulse) 
(mV/mm2)  measured  by  a  1.5mm  diameter  Titanium  probe 
versus  the  power  measured  by  a  calorimeter  per  unit  sur¬ 
face  area  (mW/mm2).  The  probe  response  is  linear  over 
a  dynamic  range  of  1.2mW  to  3.5mW  as  measured  by  the 
calorimeter.  The  error  bars  represent  3%  variation  in  pulsed 
laser  power  fluctuation  measured  by  the  probe  and  calorime¬ 
ter. 


at  0cm  from  cuvette  surface  (for  a  pathlength  in  air 
=  0.2cm).  The  beam  spread  difference  between  the 
cuvette- air  setup  v/s  cuvette- water  setup  is  accounted 
for  by  multiplying  the  fluence  measured  in  air  with  the 
ratio  of  the  beam  spreads  (or  full  width  at  half  maxi¬ 
mum  fwhm)  measured  in  the  two  settings.  The  fwhm 
values  were  measured  for  both  these  setups  using  the 
probe  (fwhmair= 40mm,  fwhmwater= 42mm)  and  their 
ratio  (0.9091)  denotes  the  proportional  change  in  fluence 
(due  to  beam-spread)  of  the  laser  beam  between  the  two 
setups.  The  experimentally  derived  calibration  factor  is 
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Fig.  6.  Filtering  the  probe  response  using  a  butterworth 
bandpass  filter.  The  left  side  shows  the  ulfiltered  signal  from 
the  probe  while  right  hand  side  shows  the  filtered  response 
with  a  high  S/N  ratio. 


then  used  for  absolute  photon  quantification  by  multi¬ 
plying  it  to  the  probe  response  measured  within  tissues 
and  tissue-like  media  (water  based  optical  phantoms)  as 
follows  (from  equation  1  and  2): 


0 absolute 


( CF  *  (^tissue)  +  1.7296 
590.46 


(3) 


Where  absolute  is  the  absolute  dose  in  mW/mm2  and 
4> tissue  is  the  fluence  measured  by  the  probe  (mV/mm2) 
in  the  tissue. 


2.D.  Optimization  of  GPU  based  Monte  Carlo 

The  original  code  of  the  GPU  based  3D  Monte  Carlo 
(MC)  was  modified  to  emit  a  broad  beam  super- 
Gaussian  beam  profile  emanating  from  an  an  optical 
fiberguide,  connected  to  a  NIR  laser.  This  profile  was 
measured  with  a  calorimeter  (gold  standard)  and  can  be 
represented  with  the  following  equation: 

4*  super  Gaussian  ^  ^ (4) 

where 

r  =  z  *  tan(6 )  (5) 

Where  z  =  perpendicular  distance  from  the  source 
to  the  plane  of  beam  profile  measurement,  0  =  equa¬ 
torial  angle  subtended  by  a  point  on  the  beam  with 
the  plane  of  the  source,  a  =  standard  deviation  of  the 
super- Gaussian  beam,  n=order  of  the  super- Gaussian 
distribution.  The  broad  beam  fiberguide  output  was 
approximated  using  multiple  single  sources  with  super- 
Gaussian  outputs.  In  our  experimental  setup,  a  match 
in  beam  profile  was  obtained  for  <7= 3.8  and  n=2.9  for 
a  fiberguide  source  of  radius  0.3cm.  Besides  the  super- 
Gaussian  beam  from  the  fiberguide,  a  flat  broad  beam 
with  radius  equal  to  the  output  port  of  the  integrat¬ 
ing  sphere  (port  radius  0.9cm)  was  used  to  characterize 
the  output  of  the  integrating  sphere.  The  beam  char¬ 
acteristics,  phantom  dimensions  and  voxel-based  optical 


properties  of  the  phantom  were  given  as  input  to  the 
MC  code  to  emulate  the  illumination  conditions  used  in 
the  optical  dosimetry  setup  of  figure  4(B). 

2.E.  Validation  of  3D  Monte  Carlo  in  brain  phan¬ 
toms  using  Optical  dosimetry 

As  described  in  section  2. A,  liquid  optical  phantoms, 
resembling  white  matter,  gray  matter  and  skull  bone, 
were  designed  using  predetermined  concentrations  of  In¬ 
dia  ink  and  Intralipid.  The  setup  for  the  light  dosimetry 
in  the  phantoms  is  similar  to  that  shown  in  figure  4(B). 
Two  sources  of  NIR  light  (an  integrating  sphere  and 
an  optical  fiberguide)  were  used  to  illuminate  the  phan¬ 
toms,  both  of  which  were  connected  to  the  laser  source 
(same  specifications  as  used  in  the  calibration  process). 
The  dosimetry  probe  was  inserted  into  the  phantom  and 
translated  using  a  calibrated  translation  stage  to  mea¬ 
sure  the  optical  fluence  (proportional  to  the  output  volt¬ 
age)  distribution  along  x,  y  and  z  directions  within  the 
phantom.  The  x,  y  and  z  fluence  distribution  generated 
by  the  Monte  Carlo  was  compared  to  the  probe  mea¬ 
surements  to  the  test  the  validity  and  accuracy  of  the 
Monte  Carlo  code. 

In  order  to  quantify  the  photon  energy  fluence  in 
J/cm 2  brain  tissues,  we  converted  the  output  of  the 
Monte  Carlo  fluence  (</>,  units  :  photons /mm2)  as  fol¬ 
lows: 


^dose{J/mm2)  =  0*—  (6) 

Where  h  is  planck’s  constant  (6.626  *  10  -34Js),  c  and 
A  are  the  speed  and  wavelength  of  light  in  the  medium. 
The  energy  deposited  in  the  tissue  is  obtained  by  mul¬ 
tiplying  the  photon  dose  with  the  absorption  coefficient 
pa  and  the  volume  of  the  voxel  (v): 


(frabsorbed^J / 'bU'bR  )  —  0  *  Pa  *  ^  ^  (7) 

2.F.  A  fast  Empirical  algorithm  for  photon  propaga¬ 
tion 

The  Empirical  approach  is  based  on  the  assumption  that 
the  photon  fluence  in  a  voxel  in  a  particular  layer  is  a 
weighted  sum  of  the  fluences  of  neighboring  voxels  in  the 
previous  layer  [17].  The  schematic  of  this  algorithm  is 
shown  in  figure  7.  The  steps  of  the  algorithm  can  be 
summarized  as  follows: 

1)  Assign  the  incident  photon  beam  distribution  (e.g. 
super- Gaussian  beam)  to  the  first  layer  of  voxels.  In  this 
case,  we  assume  that  the  direction^  of  photon  propaga¬ 
tion  is  from  layer  z=l  to  z=dim3. 

2)  For  every  voxel  in  a  layer,  calculate  the  14  direc¬ 
tional  scatter  components  from  the  incident  photon  flu- 
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ence  from  the  previous  layer: 

If(x,y,z )  =  If(x,y,z  —  l)Wtl 
x+1  y+ 1 

+  ^^{Ifd(x,y,z  —  l)Wt2) 

x-ly-1  (8) 

X  +  1  y+l 

+ X]  X+c*-^  z)wt3) 

x—l y— 1 

x+1  y+l 

Ifd(x,y,z)  =  fd(x,y,z- l)Wtl) 

x—l y— 1 
x+1  y+l 

+  X  (9) 

x—l y— 1 
x+1  y+l 

+  XX(W  z)wt2) 

x—l y— 1 

x+1  y+l 

Is(x,y,z)  =  EE(W 
x—l y— 1 
x+1  y+l 

+  ^(JMx,y„-l)^2) 

x—l y— 1 
x+1  y+l 

+  X  X+(*.^  _i)iet3)  (10) 

x—l y— 1 
x+1  y+l 

A  ^  E(V,y++D^3) 

x—l y— 1 
x+1  y+l 

+  X  ^2dbd(x,y,  z+l)Wt2) 

x—l y— 1 

Ib(x,y,z)  =  Ib(x,y,z+ l)W)tl 
x+1  y+l 

+  X  Y^dbd{x,y,z+l)Wt2) 

x-ly-1  (11) 

x+1  y+l 

+  X  X+(z>y-  z)wt3) 

x—l y— 1 

x+1  y+l 

Ibd(x,y,z)  ^  ^  ^^(-^6<i(x,y, £+1)11^1) 
x—l y— 1 
x+1  y+l 

+  X  X(+  x,y,^  +1)^2)  (12) 

x—l y— 1 
x+1  y+l 

+  XX(W  Z)W t2) 

x—l y— 1 


I(x,y,z)  ~  If(x,y,z)  +  X(-^ fd(x,y,z)) 

+  X(+^.  2))+^(x,y,^)  (18) 

"b  ^  ]( Ibd(x,y,z )) 

The  14  scatter  components  are:  1  forward  (//)  com¬ 
ponent,  1  backward  (/&)  component,  4  side  (Js)  com¬ 
ponents  and  8  diagonal  components  (4  forward  diago¬ 
nals  (Ifd)  and  4  backward  diagonals  (hd))-  The  arbi¬ 
trary  factors  wtl  and  wt2  are  the  forward  and  diagonal 
weights  assigned  to  each  directional  component,  which 
determines  the  percentage  of  scatter  and  the  sum  of  all 
weights  is  always  equal  to  one  (i.e.  wtl+4wt4=l),  which 
only  allows  for  redistribution  of  fluence  (total  fluence 
within  a  voxel  remains  same).  The  sum  of  all  the  scat¬ 
ter  components  in  a  voxel  equals  the  total  fluence  in  that 
voxel  (arbitrary  units:  photons/unit  volume). 

3)  Calculate  absorbed  fluence  (Iabs)  using  Beer  Lam¬ 
bert’s  law  for  a  pathlength  dl  =  voxel  length,  based  on 
the  absorption  coefficient  (/ia)  of  the  voxel: 

Iabs(x,y,z)  =  I(X,y,z)  *  (1  -  (14) 

4)  Calculate  and  accumulate  the  fluence  over  all  voxels 
in  the  current  layer  and  increment  z  to  proceed  to  the 
next  layer.  Continue  propagation  till  the  fluence  exits 
the  volume. 

The  GPU  based  Monte  Carlo  was  used  as  the  golden 
standard  to  derive  the  weights,  which  depend  on  the 
reduced  scattering  coefficient  of  the  medium.  The 
weights  to  be  assigned  for  scatter  were  obtained  by  an 
iterative  optimization  phase  which  compares  the  output 
fluence  of  the  empirical  algorithm  to  that  generated 
by  the  Monte  Carlo  routine  along  the  beam  axis.  The 
weights  which  show  the  best  possible  match  (using 
highest  R-square  coefficient)  are  selected  to  calibrate 
the  empirical  method.  Table  2  shows  the  weights  (wtl 
and  wt2)  obtained  for  a  pencil  beam  simulation  for 
different  tissue  types.  The  R 2  values  shows  the  degree 
of  match  between  the  Monte  Carlo  and  the  Empirical 
method  along  the  beam  axis.  The  weights  obtained 
during  optimization  depend  on  the  voxel  size  and  the 
scattering  coefficient  of  the  medium.  The  absorption 
coefficient  is  scaled  to  the  voxel  size  and  the  photons  are 
absorbed  before  proceeding  to  the  next  layer.  Figure 
11  shows  an  example  of  the  match  in  fluence  profiles 
between  the  Monte  Carlo  and  the  Empirical  approach. 

3.  Results  and  Discussion 

3. A.  Optical  probe  calibration  and  photon  quantifi¬ 
cation 

The  isotropicity  (coefficient  of  variation)  of  Nylon  and 
Titanium  based  optical  dosimetry  probes  is  shown  in 
figure  5(A).  The  titanium  based  probes  showed  a  bet¬ 
ter  isotropicity  with  an  equatorial  coefficient  of  variation 
0.0386  versus  0.1929  for  Nylon  probes  and  an  azimuthal 
coefficient  of  variation  of  0.216  for  titanium  versus  0.441 
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Probe  versus  Monte  Carlo  (White  Matter,  Gray  Matter) 
(Source:  Interating  Sphere) 
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Fig.  7.  Optical  dosimetry  in  white  and  gray  matter  phantoms.  The  photon  source  used  is  an  integrating  sphere  connected  to 
a  pulsating  laser  source.  The  graphs  show  the  exponential  decrease  of  photon  fluence  along  the  beam  axis  in  white  and  gray 
matter  mimicking  phantoms.  The  Monte  Carlo  generated  fluence  closely  matches  the  fluence  measured  by  the  probe. 


Probe  versus  Monte  Carlo  (super-Gaussia  beam) 
(source:  Optical  fiberguide) 


Distance  (cm) 


•  Probe  White 
Matter 
Phantom 

♦  Probe  Gray 
Matter 
Phantom 

□  Probe  Skull 
Bone 
Phantom 

+  MC  Skull 
23July 

XMC  WM  23 
July 

—  MC  GM  23 
July 


Fig.  8.  Optical  dosimetry  in  white  matter,  gray  matter  and  skull  bone  phantoms.  The  photon  source  used  is  an  optical 
fiberguide  (connected  to  a  pulsating  laser  source)  with  a  super-Gaussian  beam  distribution.  The  Monte  Carlo  generated  fluence 
closely  matches  the  fluence  measured  by  the  probe.  This  shows  that  the  Monte  Carlo  beam  modeling  can  be  reliably  used  as 
a  reliable  estimator  of  photon  energy  distribution  in  tissue  phantoms. 


for  Nylon.  A  titanium  based  probe  with  a  1.5mm  diam¬ 
eter  was  selected  for  further  dosimetic  studies  and  had 
an  equatorial  coefficient  of  variation  of  0.0721  and  an 
azimuthal  coefficient  of  variation  0.2959.  Figure  5(B) 
shows  the  linearity  of  response  of  the  1.5mm  Titanium 
based  probe.  The  probe  response  is  linear  over  a  dy¬ 
namic  range  of  0.01074mW/mm2  to  0.031316mW/mm2 
as  measured  by  the  calorimeter.  The  standard  error 
(Standard  Deviation/ y/n)  and  percentage  error  varia¬ 
tion  in  probe  response  was  quantified  by  measuring  mul¬ 
tiple  probe  response  in  air  as  shown  in  figure  5  over 


distances  10.1  to  18.1cm  from  the  fiberguide  surface. 
The  error  in  probe  response  is  1.695%  (standard  error 
=  0.226)  and  the  %  error  in  calorimeter  readings  was 
1.1299%  (standard  error  =  0.00495).  On  the  other  hand 
the  variation  of  the  laser  light  emitted  from  the  fiber- 
guide  over  time  was  approximately  3%.  Thus  the  error 
due  to  variation  in  pulsed  optical  power  of  the  laser  was 
the  major  source  of  error  in  our  measurements.  This 
showed  that  Titanium  probes  are  better  suited  to  opti¬ 
cal  dosimetry  studies  as  compared  to  Nylon  probes  and 
hence  were  used  in  brain  phantom  studies.  While  Nylon 
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Beam  Profile  in  Air  (MC  v/s  Calorimeter) 


White  Matter  Phantom:  Lateral  Beam  profile  (Optical  Fiberguide) 
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Gray  Matter  Phantom:  Lateral  Beam  profile  (Optical  Fiberguide) 
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Skull  Phantom:  Lateral  Beam  profile  (Optical  Fiberguide) 
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Fig.  9.  Lateral  beam  profile  measurement  using  optical  dosimetry  probe  and  validation  by  Monte  Carlo.  The  laser  source 
used  was  an  optical  fiberguide  (0.3cm  diameter)  connected  to  a  pulsed  NIR  laser  emitting  at  a  wavelength  of  800nm.  (A) 
Beam  characterization  in  air  using  calorimeter.  Optical  dosimetry  to  measure  and  validate  lateral  beam  profiles  in  phantoms 
resembling:  (B)  white  matter,  (C)  gray  matter  and  (D)  skull  bone,  i.e.,  perpendicular  to  the  beam  axis.  The  measurements  have 
been  normalized  to  the  maximum  value  measured  by  the  dosimetry  probe  in  each  of  the  respective  phantoms.  The  Monte  Carlo 
generated  fluence  closely  matches  the  fluence  measured  by  the  probe.  Thus  the  beam  spread  (i.e.  optical  power  dissipation) 
due  can  be  effectively  modeled  using  Monte  Carlo. 


probes  have  been  used  before  in  dosimetry  studies,  this 
is  the  first  time  that  titanium  probes  were  designed  with 
a  better  isotropic  response. 

In  order  to  calculate  the  absolute  photon  fluence  in  tis¬ 
sues  and  tissue-like  media,  the  probe  voltage  measured 
within  the  medium  must  be  multiplied  by  the  calibra¬ 
tion  factor  to  account  for  the  change  in  photon  loss  in 
the  probe-medium  interface  with  respect  to  the  probe- 
air  interface,  where  the  calibration  was  performed.  Since 
the  refractive  index  of  tissues  is  nearly  same  as  that  of 


water,  the  calibration  factor  obtained  for  water  can  be 
used  for  biological  media.  The  probe  calibration  factor 
was  determined  to  be  1.3182.  The  technique  presented 
in  section  2C  (using  equations  1,2,3)  provides  a  simple 
method  to  determine  the  calibration  factor  experimen¬ 
tally,  and  is  a  novel  addition  to  our  research.  It  is  im¬ 
portant  to  note  that  while  equation  2  remains  same  for 
any  probe,  the  equations  1  needs  to  experimentally  de¬ 
termined,  so  as  to  derive  equation  3  for  absolute  photon 
dose  quantification. 
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Fig.  10.  Empirical  method  for  photon  propagation:  The  fluence  at  each  voxel  in  a  layer  is  a  weighted  sum  of  the  neighboring 
fluences  in  the  previous  layer.  The  model  simulates  forward  and  diagonal  scatter  using  arbitrarily  assigned  weights,  derived 
from  the  Monte  Carlo  and  which  are  in  turn  dependent  on  the  absorption  and  scattering  coefficients  of  the  medium.  The  fluence 
entering  a  voxel  is  divided  into  14  directional  components:  1  forward,  1  backward,  4  side  components,  4  forward  diagonals  and 
4  backward  diagonals. 


Table  2.  Optimization  of  Empirical  algorithm 


Tissue 

Type 

Resolution 

(mm) 

Wtl 

Wt2 

R2 

Coefficient 

White  Matter 

0.1 

0.95364 

0.01159 

1 

Gray  Matter 

0.1 

0.940728  0.014818 

1 

Skull  Bone 

0.1 

0.953512 

0.011622 

1 

Astrocytoma 

0.1 

0.933639 

0.01659 

1 

Tumor 

White  Matter 

1 

0.372639 

0.15684 

0.999916 

Gray  Matter 

1 

0.171025  0.207244 

0.999905 

Skull  Bone 

1 

0.380072  0.154982 

0.999919 

Astrocytoma 

1 

0.21174 

0.197065 

0.999935 

Tumor 

3.B.  Optical  dosimetry  and  Monte  Carlo  validation 

The  optical  fluence  measured  by  the  dosimetry  probes 
in  white  and  gray  matter  phantoms,  was  compared  to 
the  Monte  Carlo  generated  fluence.  The  design  of  op¬ 
tical  phantoms  is  crucial  for  accurately  determining  the 
localized  distribution  of  fluence  with  the  probe  before 
validating  it  using  Monte  Carlo.  The  optical  coefficients 
not  taken  into  account  in  this  study  are  the  miniscule 
scattering  properties  of  India  ink  and  the  absorption 
properties  of  intralipid.  These  effects  while  minimal  over 


a  small  range  of  concentrations,  can  be  a  source  of  uncer¬ 
tainties  over  a  larger  range  of  concentrations  in  phantom 
design. 

Figures  7  and  8  shows  the  comparison  between  the 
probe  measurements  and  the  Monte  Carlo  in  brain  phan¬ 
toms  for  two  different  illumination  profiles,  flat  diffuse 
beam  (integrating  sphere)  and  super- Gaussian  (optical 
fiberguide)  respectively.  In  addition  to  the  match  in  the 
exponential  direction,  a  comparison  was  also  made  in 
the  direction  perpendicular  to  the  beam  axis,  as  shown 
in  figure  9.  The  plots  show  that  the  fluence  measured 
by  the  probe  match  the  Monte  Carlo  simulations  both 
along  and  perpendicular  to  the  beam  direction.  The 
probe  measurements  thus  validate  the  3D  Monte  Carlo 
in  optical  phantoms  resembling  brain  tissues,  both  along 
and  perpendicular  to  the  beam  axis,  and  is  the  major 
purpose  of  this  research. 

An  analysis  of  the  dose  profile  within  the  skull  phan¬ 
tom  showed  that  the  super-Gaussian  beam  fluence  de¬ 
creases  with  an  effective  attenuation  of  4.818.  This  fac¬ 
tor  was  used  to  estimate  the  effectiveness  of  therapy  for 
different  skull  thicknesses.  Since  the  input  photon  flu¬ 
ence  of  the  fiberguide  is  approximately  13.25mW  (mea¬ 
sured  by  calorimeter  with  aperture  12mm),  the  effec¬ 
tive  photon  dose  deposited  in  a  skull  of  thickness  1cm 
is  0.107mW  (0.0072mW/rara2),  i.e.  0.808%  of  optical 
power  reaches  the  soft  tissues.  On  the  other  hand  a 
pediatric  skull  of  0.5cm  thickness  would  allow  8.99%  op¬ 
tical  power  to  reach  the  soft  tissue  (i.e.  0.0795mW/mm2 
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Fig.  11.  Comparison  of  the  Empirical  algorithm  with  the  Monte  Carlo  in  homogeneos  (a,b)  and  heterogeneous  (c,d)  brain 
phantoms  resembling  white  matter,  gray  matter  and  astrocytoma.  Part  (a)  shows  the  optimization  phase  using  pencil  beam 
simulations  while  part  2  shows  exponential  decrease  in  fluence  for  a  Gaussian  beam  simulation.  Part  (c)  shows  a  simulated 
heterogeneous  brain  phantom  while  part  (d)  shows  fluence  profiles  along  beam  axis  for  the  Empirical  and  Monte  Carlo  models. 
The  empirical  approach  is  a  close  approximation  of  the  Monte  Carlo  over  the  range  of  optical  properties  of  brain  tissues  [17]. 


from  the  fiberguide  source).  The  thickness  of  the  skull 
and  skull  density  determines  the  overall  attenuation  and 
hence  an  optical  planning  protocol  needs  an  accurate 
bone  density  and  thickness  map  which  can  be  provided 
using  a  CT  scan. 

3.C.  Empirical  versus  Monte  Carlo 

Figure  11  shows  the  comparison  between  the  Monte 
Carlo  and  Empirical  in  brain  tissues  for  both  pencil 
and  broad  beam  simulations  in  homogeneous  and  het¬ 
erogeneous  phantoms.  Figure  11(a)  shows  the  results 
for  pencil  beam  simulations,  11(b)  shows  comparison  us¬ 
ing  Gaussian  beams  while  11(d)  shows  the  comparison 
in  a  heterogeneous  simulated  head  phantom  as  shown 
in  11(c),  with  the  optical  properties  described  in  table 
1.  The  weights  derived  from  the  optimization  phase  for 
pencil  beams  were  used  in  the  broad  beam  simulations 
(figures  11(c)  and  11(d)).  The  fluence  generated  by  the 
Monte  Carlo  closely  matches  the  empirical  fluence  in 
both  homogeneous  and  heterogeneous  brain  tissue  phan¬ 
toms  along  the  beam  axis. 

The  Monte  Carlo  code  is  computationally  expensive 
in  spite  of  being  based  on  a  GPU  based  platform.  In 
order  to  simulate  the  optical  properties  and  fine  struc¬ 
tural  boundaries  of  brain  tissues,  we  need  a  fine  resolu¬ 


tion  equal  to  or  less  than  1mm.  An  optical  map  of  the 
brain  can  be  obtained  through  segmenting  and  translat¬ 
ing  CT  and/or  MR  images  into  corresponding  optical 
properties.  With  such  tiny  resolutions,  the  CT  dataset 
size  of  the  brain  becomes  very  large  (for  example  a  CT 
dataset  of  size  512x512x512  matrix  with  0.625mm  res¬ 
olution  and  160mm  collimation) .  This  significantly  in¬ 
creases  the  execution  time  of  the  Monte  Carlo  routine  to 
iterate  over  various  illumination  conditions  and  to  simu¬ 
late  photon  propagation  through  a  large  dataset  in  every 
iteration.  The  Empirical  approach  provides  an  alter¬ 
native  means  to  replace  the  computationally  expensive 
Monte  Carlo,  for  preliminary  analysis  by  reducing  the 
execution  time  to  select  the  best  illumination  conditions 
for  therapy.  The  acceleration  achieved  by  the  Empiri¬ 
cal  algorithm  was  approximately  700X  for  a  pencil  beam 
simulation  and  16000X  for  a  broad  beam  simulation  in 
a  head  phantom. 

4.  Conclusions 

This  study  validates  the  3D  fluence  distribution  simu¬ 
lated  by  a  Monte  Carlo,  using  an  optical  dosimetry  probe 
in  phantoms  resembling  white  matter,  gray  matter  and 
skull  bone.  The  use  of  a  dosimetry  probe  with  a  better 
design,  i.e.,  with  a  linear  and  highly  isotropic  response, 
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allowed  for  accurate  localized  measurement  of  NIR  pho¬ 
ton  fluence  in  tissue  phantoms  (<  1%  error),  and  fur¬ 
ther  increased  the  accuracy  of  the  measurements.  The 
probes  themselves  were  calibrated  using  calorimeter  gold 
standard  and  a  unique  calibration  factor  was  derived 
specific  to  the  probe.  The  3D  validation  of  the  Monte 
Carlo  generated  fluence  for  two  different  types  of  NIR 
light  sources  in  different  optical  phantoms,  and  the  ease 
of  photon  quantification  to  absolute  dose  (in  mW /mm2 
or  J /mm2)  proves  the  fidelity  of  the  Monte  Carlo  as  a 
predictive  tool  to  accurately  estimate  3D  photon  energy 
distribution  in  complex  media.  The  prediction  error  of 
the  Monte  Carlo  is  determined  by  the  3%  variation  of  op¬ 
tical  power  and  <  1%  error  of  the  calorimeter/dosimetry 
probe  setup.  The  voxelated  structure  of  the  medium  to 
be  simulated  allows  the  Monte  Carlo  to  emulate  com¬ 
plex  geometrical  shaped  of  tissues  by  accepting  datasets 
generated  by  other  imaging  systems  (e.g.  CT/MRI)  and 
without  any  pre-knowledge  of  tissue  geometries,  unlike 
previous  versions  of  the  Monte  Carlo.  This,  along  with 
the  flexibility  of  generating  user  defined  photon  beams 
proves  adaptability  of  the  Monte  Carlo  to  be  used  in 
clinical  studies  as  a  predictive  tool  to  estimate  photon 
dose  distribution  in  heterogeneous  brain  tissues. 

The  empirical  approach  provides  a  robust  and  fast 
means  of  assessing  the  tissue  illumination  and  its  weights 
are  obtained  through  the  optimization  phase  with  the 
Monte  Carlo.  The  acceleration  provided  by  the  Empir¬ 
ical  algorithm  increases  the  clinical  feasibility  of  optical 
therapeutic  planning  to  narrow  down  the  complex  pos¬ 
sibilities  of  illumination  conditions,  further  compounded 
by  the  heterogeneous  structure  of  the  brain  (e.g.  varying 
skull  thicknesses  and  densities).  One  disadvantage  of  us¬ 
ing  the  Empirical  approach  is  that  the  weights  assigned 
depend  on  the  resolution  of  the  medium  and  changing 
resolution  requires  an  optimization  phase  to  derive  the 
weights.  In  spite  of  its  limitations  and  optimization  re¬ 
quirements,  the  Empirical  can  provide  an  acceleration 
over  the  Monte  Carlo  to  analyze  complex  illuminations 
via  iterative  routines  and  in  preliminary  studies  to  de¬ 
termine  the  best  irradiation  angles.  The  empirical  ap¬ 
proach  can  thus  be  used  as  a  preliminary  step  prior  to 
using  the  Monte  Carlo,  with  the  final  analysis  (final  es¬ 
timation)  being  done  by  the  Monte  Carlo  to  ensure  high 
therapeutic  accuracy  to  ensure  drug  release. 

Our  ultimate  goal  is  to  design  a  fast  Monte  Carlo 
based  optical  therapeutic  protocol  to  treat  brain  metas¬ 
tasis.  This  requires  voxel  by  voxel  determination  of 
photon  energy  distribution  to  accurately  predict  the 
rate/quantity  of  drug  release.  This  study  achieves  the 
first  step  of  designing  such  a  tool  to  validate  the  NIR 
photon  distribution  in  three  dimensions  in  optical  phan¬ 
toms.  The  next  step  will  be  to  test  Monte  Carlo  in 
predicting  photon  distribution  in  heterogeneous  tissue 
samples  and  under  various  illuminations  beam  profiles. 
The  tissue  structure  and  boundaries  can  be  determined 
using  other  imaging  modalities  such  as  MRI  or  CT  scan 
to  separate  various  tissue  types  and  obtain  an  optical 


property  map  of  the  tissue  sample.  This  can  be  used  an 
input  to  the  Monte  Carlo  to  obtain  photon  distribution 
in  the  tissue  sample. 
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